Seberapa berarti, atau berguna, kesalahan parameter yang dihasilkan saat melakukan LinearModelFit atau NonlinearModelFit yang tidak berbobot?

Aug 16 2020

Ini mungkin pertanyaan yang tergelincir di ambang milik lebih ke bidang statistik dan SE Cross-Validated , tapi saya juga secara khusus tertarik dengan rutinitas pemasangan Mathematica.

Biasanya, jika saya ingin menyesuaikan model ke beberapa data menggunakan salah satunya NonlinearModelFitatau LinearModelFitsaya akan mengalami beberapa kesalahan yang terkait dengan file$y$-data yang saya gunakan untuk menimbang pas. Kesalahan ini mungkin hanya kesalahan standar yang diperoleh dari pengukuran berulang, atau saya mungkin tahu sesuatu tentang proses fisik dan dapat menetapkan bobot.

Misalnya Weights->1/YDataErrors^2dan saya selalu menetapkan penaksir varians saya sebagai VarianceEstimatorFunction -> (1 &). Saya kemudian bisa mendapatkan kesalahan parameter saya dari matriks kovarians, atau cukup dengan MyFit["ParameterErrors"].

Namun dalam beberapa kasus seseorang mungkin tidak memiliki kesalahan untuk data yang ingin Anda paskan, yang berarti seseorang tidak dapat memberikan bobot seperti yang saya jelaskan di atas. Pertanyaan saya kemudian, seberapa dapat diandalkan - atau yang lebih penting - seberapa bermakna secara fisik / statistik adalah kesalahan parameter untuk kesesuaian tak berbobot di Mathematica?

Jawaban

4 JimB Aug 16 2020 at 19:20

Misalnya, jika seseorang memiliki dua sumber kesalahan, katakan kesalahan pengukuran dan kesalahan kurangnya kesesuaian, maka menggunakan bobot berdasarkan kesalahan pengukuran dapat menghasilkan perkiraan yang terlalu rendah dari kesalahan standar. Pertimbangkan model berikut:

$$y=a+b x +\gamma + \epsilon$$

dimana $y$ adalah respon yang diukur, $x$ adalah prediktornya, $a$ dan $b$ adalah konstanta untuk diperkirakan, $\gamma$ adalah kesalahan pengukuran berulang dengan $\gamma \sim N(0,\sigma_{ME})$, dan $\epsilon$ adalah kesalahan kurang pas dengan $\epsilon \sim N(0,\sigma)$ dan semua kesalahan dianggap independen.

Pertama, setel beberapa parameter tertentu:

(* Measurement error standard deviation *)
σME = 10;

(* Lack-of-fit error standard deviation *)
σ = 20;

(* Regression coefficients *)
a = 1;
b = 1;

Hasilkan dan plot beberapa data:

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]

Sekarang pertimbangkan dua model linier yang berbeda cocok di mana lm1apa yang Anda sarankan dan lm2apa yang saya sarankan:

lm1 = LinearModelFit[data, z, z, Weights -> 1/ConstantArray[σME^2, n],
   VarianceEstimatorFunction -> (1 &)];
lm2 = LinearModelFit[data, z, z];
lm1["ParameterTable"]

lm2["ParameterTable"]

Estimasi parameternya identik tetapi error standarnya lm1kurang dari setengah ukurannya lm2. Yang mana yang benar?

Matriks kovariansi yang "benar" dari penduga kuadrat terkecil dari adan buntuk model ini adalah

$$\left(\sigma ^2+\sigma_{ME}^2\right) \left(X^T.X\right)^{-1}$$

dimana $X$adalah matriks desain. Dalam kode Mathematica kesalahan standarnya badalah

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

Itu sangat cocok dengan lm2.

Ini adalah contoh yang sedikit dibuat-buat karena saya memiliki semua kesalahan standar pengukuran yang identik karena fungsi regresi Mathematica hanya mengizinkan satu istilah kesalahan. Dan dengan memiliki kesalahan standar pengukuran yang identik, yang menghasilkan model yang setara dengan kesalahan tunggal.

Namun, meskipun deviasi standar pengukuran sangat bervariasi, masalah tentang pembobotan yang tidak tepat sehingga tidak sesuai dengan struktur kesalahan model tetap ada.

Rutinitas regresi Mathematica belum memadai untuk model dengan lebih dari satu sumber kesalahan. Saya berharap mereka begitu.