Насколько значимы или полезны ошибки параметров, возникающие при выполнении невзвешенных LinearModelFit или NonlinearModelFit?

Aug 16 2020

Это может быть вопрос, который балансирует на грани того, чтобы больше относиться к сфере статистики и перекрестной проверки SE , но меня также особенно интересуют процедуры подгонки Mathematica.

Обычно, если я хочу подогнать модель к некоторым данным, используя либо, NonlinearModelFitлибо, у LinearModelFitменя будут ошибки, связанные с моим$y$-данные, которые я использую для взвешивания припадков. Эти ошибки могут быть просто стандартной ошибкой, полученной при повторных измерениях, или я могу кое-что знать о физических процессах и могу назначать веса.

Например, Weights->1/YDataErrors^2я всегда устанавливаю свою оценку дисперсии как VarianceEstimatorFunction -> (1 &). Затем я могу получить ошибки параметров из ковариационной матрицы или просто с помощью MyFit["ParameterErrors"].

Однако в некоторых случаях у вас может не быть ошибок для данных, которые вы хотите подогнать, что означает, что нельзя предоставить веса так, как я описал выше. Тогда мой вопрос: насколько надежны - или, что более важно - насколько физически / статистически значимы ошибки параметров для невзвешенной подгонки в системе Mathematica?

Ответы

4 JimB Aug 16 2020 at 19:20

Например, если имеется два источника ошибок, скажем, ошибка измерения и ошибка несовпадения, то использование весов, основанных на ошибках измерения, может привести к существенному занижению стандартных ошибок. Рассмотрим следующую модель:

$$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 коде стандартная ошибка bIS

X = Transpose[{ConstantArray[1, n], Range[n]}]
Sqrt[(σME^2 + σ^2) Inverse[Transpose[X].X][[2, 2]]] // N
(* 0.0774635 *)

Это очень хорошо сочетается с lm2.

Это слегка надуманный пример, поскольку у меня все стандартные ошибки измерений идентичны, потому что функции регрессии Mathematica допускают только один член ошибки. И при идентичности стандартных ошибок измерения получается эквивалентная модель с единственной ошибкой.

Однако даже когда стандартные отклонения измерений значительно различаются, остается проблема неправильного взвешивания, которое не соответствует структуре ошибок модели.

Подпрограммы регрессии в системе Mathematica еще не подходят для моделей с более чем одним источником ошибок. Я бы хотел, чтобы они были.