重みによるモデル近似の改善
この例では、比較のために線形最小二乗法と重み付き最小二乗法の両方を使用して多項式モデルをデータに当てはめる方法を説明します。
関数randn
を使用して、さまざまな正規分布からサンプル データを生成します。
rng("default") % For reproducibility rnorm = []; idx = []; for k=1:20 r = k*randn([20,1]) + (1/20)*(k^3); rnorm = [rnorm;r]; idx = [idx;ones(20,1).*k]; end
従属変数 rnorm
には、20 個の正規分布からのサンプル データが含まれます。独立変数 idx
には、rnorm
内の 2 つの要素が同じ正規分布からサンプリングされているかどうかを示す整数が含まれます。
3 次多項式モデルを idx
と rnorm
に当てはめます。モデルを当てはめるときに使用された係数推定値とアルゴリズムに関する情報を返します。
[p3fit,~,fitinfo] = fit(idx,rnorm,"poly3")
p3fit = Linear model Poly3: p3fit(x) = p1*x^3 + p2*x^2 + p3*x + p4 Coefficients (with 95% confidence bounds): p1 = 0.05156 (0.0438, 0.05932) p2 = -0.03993 (-0.2875, 0.2076) p3 = 0.1418 (-2.124, 2.408) p4 = 0.0462 (-5.585, 5.678)
fitinfo = struct with fields:
numobs: 400
numparam: 4
residuals: [400x1 double]
Jacobian: [400x4 double]
exitflag: 1
algorithm: 'QR factorization and solve'
iterations: 1
p3fit
には、モデル係数 (95% 信頼区間) の推定値が含まれます。多項式モデルを当てはめる場合の既定の近似法は線形最小二乗法です。fitinfo
には、係数をデータに当てはめるときに使用された近似アルゴリズムについての情報が含まれます。データ内の誤差は、fitinfo
に保存される残差によって推定できます。
ステム プロットを使用して残差をプロットします。
stem(idx,fitinfo.residuals) xlabel("idx") ylabel("residuals")
残差のプロットは、idx
内の値が大きくなるにつれて、分散が大きくなっていることを示しています。idx
のさまざまな値に対して分散が一定でないことは、重み付き最小二乗近似法の方がモデル係数の計算に適していることを示しています。
後で重みを保存するために、ゼロのベクトルを作成します。
W = zeros(length(rnorm),1);
ここで提供する重みは、残差分散を変換し、それが idx
のさまざまな値に対して一定になるようにします。rnorm
内の要素ごとに、idx
内の対応する値の残差分散の逆数として重みを定義します。次に、重み付きのモデルを当てはめます。
for k=1:20 rnorm_idx = rnorm(idx==k); recipvar = 1/var(rnorm_idx); w = (idx==k).*recipvar; W = W+w; end [wp3fit,~,wfitinfo] = fit(idx,rnorm,"poly3","Weights",W)
wp3fit = Linear model Poly3: wp3fit(x) = p1*x^3 + p2*x^2 + p3*x + p4 Coefficients (with 95% confidence bounds): p1 = 0.04894 (0.04419, 0.0537) p2 = 0.03601 (-0.08744, 0.1595) p3 = -0.4262 (-1.253, 0.4009) p4 = 0.9836 (-0.1959, 2.163)
wfitinfo = struct with fields:
numobs: 400
numparam: 4
residuals: [400x1 double]
Jacobian: [400x4 double]
exitflag: 1
algorithm: 'QR factorization and solve'
iterations: 1
wp3fit
と wfitinfo
には、重み付き最小二乗近似の結果が含まれます。
p3fit
、wp3fit
、rnorm
を同じプロットに表示します。
plot(p3fit,idx,rnorm) hold on plot(wp3fit,"g") xlabel("idx") ylabel("rnorm") legend(["rnorm","linear least-squares fit","weighted least-squares fit"]) hold off
プロットは、wp3fit
が p3fit
を密接に追従していることを示しています。
wp3fit
の方が p3fit
より適切な近似かどうかは、残差をプロットすることによって判断できます。
stem(idx,wfitinfo.residuals) xlabel("idx") ylabel("residuals")
この出力は、wp3fit
残差の方が、p3fit
残差より小さいことを示しています。また、wp3fit
残差の分散の方が、p3fit
残差の分散より、idx
の値が変わっても類似しています。