Hmm. Not sure how this is relevant to a linear regression?

No, the variance in each bin is not a new predictor, but would combine with the variance in the dependent variables to give a weight in the chisquare sum (in the way which I indicated above, maybe there should be a square root). This seems intuitively true, but not sure how rigorous it is as far as obtaining true parameters from the limit of infinite data sets. Is there a Mathematica function or option to LinearModelFit that can do this?

I couldn't find these options with an online search here for LinearModelFit

just some comment in defense of the mathematicas:

1. stepping stones in augmenting power of statistics functions were V8 and V10. statistics functionality was heavily revised there, so you're missing out on the augmentation of V10.
2. by now there are even more statistics functions, see the Repository launched in summer 2019.
3. it may or may not be true that $R$ or $Maple$ have even more of those specialized statistics functions but there is some chance that none has the exact special function that you're looking for. If you do find the exact function there which you are looking for, then please share a link and we will understand exactly what you're looking for. You're request is too vague, i doht understand it completely without an example, sorry.
4. $\text{Wolfram L}$ provides you all the functions and vocab items you need to create your own specialized statistics functions (utility functions), tailored to your needs. So if you need to compose your own function (because none of the competition has it), then you're in the right programming environment. And if statistics is all you do, day in day out, then a specialized comprehensive GUI-based statistics application would be the better choice.

Keep up posted how your quest is going!

Hi!

Notice that "TotalVariation" and "Tikhonv" are options for Fit. But I wonder if there is a way to determine the regularization parameter?

https://reference.wolfram.com/language/ref/ParameterMixtureDistribution.html

I'm going to try to revive this post. If there are known variances in dependent variable, their inverse can be used as a weight in computing the chi square sum in the regression analysis. But what if we also have information about the independent variables' variances? Can we use that information somehow? Naively, I would say (for one independent variable xi and one dependent variable yi) to find the regression coefficient r that minimizes chisquare=Sumi[(r xi -yi)^2/(r^2 sigmaxi^2+sigmay_i^2)]. Is this correct? And if so, how can the Mathematica stat package implement it (I suppose I could write my own code to do it, if not available in stat package)