Насколько значимы или полезны ошибки параметров, возникающие при выполнении невзвешенных LinearModelFit или NonlinearModelFit?
Это может быть вопрос, который балансирует на грани того, чтобы больше относиться к сфере статистики и перекрестной проверки SE , но меня также особенно интересуют процедуры подгонки Mathematica.
Обычно, если я хочу подогнать модель к некоторым данным, используя либо, NonlinearModelFit
либо, у LinearModelFit
меня будут ошибки, связанные с моим$y$-данные, которые я использую для взвешивания припадков. Эти ошибки могут быть просто стандартной ошибкой, полученной при повторных измерениях, или я могу кое-что знать о физических процессах и могу назначать веса.
Например, Weights->1/YDataErrors^2
я всегда устанавливаю свою оценку дисперсии как VarianceEstimatorFunction -> (1 &)
. Затем я могу получить ошибки параметров из ковариационной матрицы или просто с помощью MyFit["ParameterErrors"]
.
Однако в некоторых случаях у вас может не быть ошибок для данных, которые вы хотите подогнать, что означает, что нельзя предоставить веса так, как я описал выше. Тогда мой вопрос: насколько надежны - или, что более важно - насколько физически / статистически значимы ошибки параметров для невзвешенной подгонки в системе Mathematica?
Ответы
Например, если имеется два источника ошибок, скажем, ошибка измерения и ошибка несовпадения, то использование весов, основанных на ошибках измерения, может привести к существенному занижению стандартных ошибок. Рассмотрим следующую модель:
$$y=a+b x +\gamma + \epsilon$$
где $y$ это измеренный отклик, $x$ предсказатель, $a$ и $b$ константы, которые необходимо оценить, $\gamma$ повторная ошибка измерения с $\gamma \sim N(0,\sigma_{ME})$, и $\epsilon$ ошибка несоответствия с $\epsilon \sim N(0,\sigma)$ и все ошибки считаются независимыми.
Сначала установите некоторые конкретные параметры:
(* Measurement error standard deviation *)
σME = 10;
(* Lack-of-fit error standard deviation *)
σ = 20;
(* Regression coefficients *)
a = 1;
b = 1;
Сгенерируйте и нанесите на график некоторые данные:
n = 100;
x = Range[n];
SeedRandom[12345];
measurementError = RandomVariate[NormalDistribution[0, σME], n];
lackOfFitError = RandomVariate[NormalDistribution[0, σ], n];
y = a + b x + measurementError + lackOfFitError;
data = Transpose[{x, y}];
data2 = {#[[1]], Around[#[[2]], σME]} & /@ data;
ListPlot[data2]
Теперь рассмотрим две разные линейные модели, где lm1
это то, что вы предлагаете, и lm2
то, что предлагаю я:
lm1 = LinearModelFit[data, z, z, Weights -> 1/ConstantArray[σME^2, n],
VarianceEstimatorFunction -> (1 &)];
lm2 = LinearModelFit[data, z, z];
lm1["ParameterTable"]
lm2["ParameterTable"]
Оценки параметров идентичны, но стандартные ошибки для lm1
меньше половины размера, чем для lm2
. Который правильный?
«Истинная» ковариационная матрица оценок наименьших квадратов для этой модели a
и b
для этой модели равна
$$\left(\sigma ^2+\sigma_{ME}^2\right) \left(X^T.X\right)^{-1}$$
где $X$это матрица дизайна. В Mathematica коде стандартная ошибка b
IS
X = Transpose[{ConstantArray[1, n], Range[n]}]
Sqrt[(σME^2 + σ^2) Inverse[Transpose[X].X][[2, 2]]] // N
(* 0.0774635 *)
Это очень хорошо сочетается с lm2
.
Это слегка надуманный пример, поскольку у меня все стандартные ошибки измерений идентичны, потому что функции регрессии Mathematica допускают только один член ошибки. И при идентичности стандартных ошибок измерения получается эквивалентная модель с единственной ошибкой.
Однако даже когда стандартные отклонения измерений значительно различаются, остается проблема неправильного взвешивания, которое не соответствует структуре ошибок модели.
Подпрограммы регрессии в системе Mathematica еще не подходят для моделей с более чем одним источником ошибок. Я бы хотел, чтобы они были.