Teori di balik argumen bobot dalam R saat menggunakan lm ()

Setelah satu tahun di sekolah pascasarjana, pemahaman saya tentang "kuadrat terkecil tertimbang" adalah sebagai berikut: biarkan $\mathbf{y} \in \mathbb{R}^n$ , $\mathbf{X}$ menjadi beberapa matriks desain $n \times p$ , menjadi parameter vektor, \ boldsymbol \ epsilon \ in \ mathbb {R} ^ n menjadi vektor kesalahan sedemikian rupa sehingga \ boldsymbol \ epsilon \ sim \ mathcal {N} (\ mathbf {0}, \ sigma ^ 2 \ mathbf {V}) , di mana \ mathbf {V} = \ text {diag} (v_1, v_2, \ dots, v_n) dan \ sigma ^ 2> 0 . Kemudian model \ mathbf {y} = \ mathbf {X} \ boldsymbol \ beta + \ boldsymbol \ epsilon$\boldsymbol\beta \in \mathbb{R}^p$$\boldsymbol\epsilon \in \mathbb{R}^n$$\boldsymbol\epsilon \sim \mathcal{N}(\mathbf{0}, \sigma^2\mathbf{V})$$\mathbf{V} = \text{diag}(v_1, v_2, \dots, v_n)$$\sigma^2 > 0$

$$\mathbf{y} = \mathbf{X}\boldsymbol\beta + \boldsymbol\epsilon$$ di bawah asumsi disebut model "kuadrat terkecil". Masalah WLS pada akhirnya adalah menemukan Misalkan \ mathbf {y} = \ begin {bmatrix } y_1 & \ dots & y_n \ end {bmatrix} ^ {T} , \ boldsymbol \ beta = \ begin {bmatrix} \ beta_1 & \ dots & \ beta_p \ end {bmatrix} ^ {T} , dan \ mathbf {X } = \ begin {bmatrix} x_ {11} & \ cdots & x_ {1p} \\ x_ {21} & \ cdots & x_ {2p} \\ \ vdots & \ vdots & \ vdots \\ x_ {n1} & & \ cdots & x_ {np} \ end {bmatrix} = \ begin {bmatrix} \ mathbf {x} _ {1} ^ {T} \\ \ mathbf {x} _ {2} ^ {T} \\ \ vdots \\ \ mathbf {x} _ {n} ^ {T} \ end {bmatrix} \ text {. $$\begin{equation} \arg\min_{\boldsymbol \beta}\left(\mathbf{y}-\mathbf{X}\boldsymbol\beta\right)^{T}\mathbf{V}^{-1}\left(\mathbf{y}-\mathbf{X}\boldsymbol\beta\right)\text{.} \end{equation}$$ $\mathbf{y} = \begin{bmatrix} y_1 & \dots & y_n\end{bmatrix}^{T}$$\boldsymbol\beta = \begin{bmatrix} \beta_1 & \dots & \beta_p\end{bmatrix}^{T}$ $$\mathbf{X} = \begin{bmatrix} x_{11} & \cdots & x_{1p} \\ x_{21} & \cdots & x_{2p} \\ \vdots & \vdots & \vdots \\ x_{n1} & \cdots & x_{np} \end{bmatrix} = \begin{bmatrix} \mathbf{x}_{1}^{T} \\ \mathbf{x}_{2}^{T} \\ \vdots \\ \mathbf{x}_{n}^{T} \end{bmatrix}\text{.}$$ $\mathbf{x}_i^{T}\boldsymbol\beta\in \mathbb{R}^1$ , jadi $$\mathbf{y}-\mathbf{X}\boldsymbol\beta = \begin{bmatrix} y_1-\mathbf{x}_{1}^{T}\boldsymbol\beta \\ y_2-\mathbf{x}_{2}^{T}\boldsymbol\beta \\ \vdots \\ y_n-\mathbf{x}_{n}^{T}\boldsymbol\beta \end{bmatrix}\text{.}$$ Ini memberi $$\begin{align} (\mathbf{y}-\mathbf{X}\boldsymbol\beta)^{T}\mathbf{V}^{-1} &= \begin{bmatrix} y_1-\mathbf{x}_{1}^{T}\boldsymbol\beta &y_2-\mathbf{x}_{2}^{T}\boldsymbol\beta & \cdots & y_n-\mathbf{x}_{n}^{T}\boldsymbol\beta \end{bmatrix}\text{diag}(v_1^{-1}, v_2^{-1}, \dots, v_n^{-1}) \\ &= \begin{bmatrix} v_1^{-1}(y_1-\mathbf{x}_{1}^{T}\boldsymbol\beta) &v_2^{-1}(y_2-\mathbf{x}_{2}^{T}\boldsymbol\beta) & \cdots & v_n^{-1}(y_n-\mathbf{x}_{n}^{T}\boldsymbol\beta) \end{bmatrix} \end{align}$$ v_n ^ {- 1} (y_n- \ mathbf {x} _ {n} ^ {T} \ boldsymbol \ beta) \ end {bmatrix} \ end {align} sehingga memberikan $$\begin{equation} \arg\min_{\boldsymbol \beta}\left(\mathbf{y}-\mathbf{X}\boldsymbol\beta\right)^{T}\mathbf{V}^{-1}\left(\mathbf{y}-\mathbf{X}\boldsymbol\beta\right) \end{equation} = \arg\min_{\boldsymbol \beta}\sum_{i=1}^{n}v_i^{-1}(y_i-\mathbf{x}^{T}_i\boldsymbol\beta)^2\text{.}$$ $\boldsymbol\beta$ diperkirakan menggunakan $$\hat{\boldsymbol\beta} = (\mathbf{X}^{T}\mathbf{V}^{-1}\mathbf{X})^{-1}\mathbf{X}^{T}\mathbf{V}^{-1}\mathbf{y}\text{.}$$ Ini adalah tingkat pengetahuan yang saya kenal. Saya tidak pernah diajari bagaimana $v_1, v_2, \dots, v_n$ harus dipilih, walaupun tampaknya, dilihat dari sini , biasanya $\text{Var}(\boldsymbol\epsilon) = \text{diag}(\sigma^2_1, \sigma^2_2, \dots, \sigma^2_n)$, yang masuk akal secara intuitif. (Berikan bobot yang sangat bervariasi, lebih sedikit bobot dalam masalah WLS, dan berikan pengamatan dengan variabilitas yang lebih sedikit, lebih berat.)

Apa yang saya sangat ingin tahu tentang bagaimana Rmenangani bobot dalam lm()fungsi ketika bobot ditugaskan menjadi bilangan bulat. Dari menggunakan ?lm:

Non- NULLbobot dapat digunakan untuk menunjukkan bahwa pengamatan yang berbeda memiliki varian yang berbeda (dengan nilai dalam bobot berbanding terbalik dengan varian); atau ekuivalen, ketika elemen bobot adalah bilangan bulat positif , bahwa setiap respons adalah rata-rata dari satuan berat (termasuk kasus bahwa ada pengamatan sama dengan dan data telah dirangkum).$w_i$$y_i$$w_i$$w_i$$y_i$

Saya telah membaca ulang paragraf ini beberapa kali, dan itu tidak masuk akal bagi saya. Menggunakan kerangka kerja yang saya kembangkan di atas, misalkan saya memiliki nilai simulasi berikut:

x <- c(0, 1, 2)
y <- c(0.25, 0.75, 0.85)
weights <- c(50, 85, 75)

lm(y~x, weights = weights)

Call:
lm(formula = y ~ x, weights = weights)

Coefficients:
(Intercept)            x  
     0.3495       0.2834  

Menggunakan kerangka yang saya kembangkan di atas, bagaimana parameter ini diturunkan? Inilah upaya saya untuk melakukan ini dengan tangan: dengan asumsi , kami telah dan melakukan ini dalam memberi (perhatikan bahwa keterbalikan tidak berfungsi dalam kasus ini, jadi saya menggunakan invers umum):$\mathbf{V} = \text{diag}(50, 85, 75)$

$$\begin{align}&\begin{bmatrix} \hat\beta_0 \\ \hat\beta_1 \end{bmatrix} = \\ &\left(\begin{bmatrix} 1 & 1\\ 1 & 1\\ 1 & 1 \end{bmatrix}\text{diag}(1/50, 1/85, 1/75)\begin{bmatrix} 1 & 1\\ 1 & 1\\ 1 & 1 \end{bmatrix}^{T} \right)^{-1}\begin{bmatrix} 1 & 1\\ 1 & 1\\ 1 & 1 \end{bmatrix}^{T}\text{diag}(1/50, 1/85, 1/75)\begin{bmatrix} 0.25 \\ 0.75 \\ 0.85 \end{bmatrix} \end{align}$$ R

X <- matrix(rep(1, times = 6), byrow = T, nrow = 3, ncol = 2)
V_inv <- diag(c(1/50, 1/85, 1/75))
y <- c(0.25, 0.75, 0.85)

library(MASS)
ginv(t(X) %*% V_inv %*% X) %*% t(X) %*% V_inv %*% y

         [,1]
[1,] 0.278913
[2,] 0.278913

Ini tidak cocok dengan nilai dari lm()output. Apa yang saya lakukan salah?

Matriks harus bukan Anda juga harus , bukan .$X$

$$\begin{bmatrix} 1 & 0\\ 1 & 1\\ 1 & 2 \end{bmatrix},$$ $$\begin{bmatrix} 1 & 1\\ 1 & 1\\ 1 & 1 \end{bmatrix}.$$ V_invdiag(weights)diag(1/weights)

x <- c(0, 1, 2)
y <- c(0.25, 0.75, 0.85)
weights <- c(50, 85, 75)
X <- cbind(1, x)

> solve(t(X) %*% diag(weights) %*% X, t(X) %*% diag(weights) %*% y)
       [,1]
  0.3495122
x 0.2834146

Untuk menjawab ini dengan lebih ringkas, regresi kuadrat terkecil yang digunakan weightsdalam Rmembuat asumsi berikut: misalkan kita miliki weights = c(w_1, w_2, ..., w_n). Biarkan , menjadi matriks desain, menjadi vektor parameter, dan menjadi vektor kesalahan dengan mean dan matriks varian , di mana . Kemudian, Mengikuti langkah-langkah derivasi yang sama dalam posting asli, kami memiliki $\mathbf{y} \in \mathbb{R}^n$$\mathbf{X}$$n \times p$$\boldsymbol\beta\in\mathbb{R}^p$$\boldsymbol\epsilon \in \mathbb{R}^n$$\mathbf{0}$$\sigma^2\mathbf{V}$$\sigma^2 > 0$

$$\mathbf{V} = \text{diag}(1/w_1, 1/w_2, \dots, 1/w_n)\text{.}$$ $$\begin{align} \arg\min_{\boldsymbol \beta}\left(\mathbf{y}-\mathbf{X}\boldsymbol\beta\right)^{T}\mathbf{V}^{-1}\left(\mathbf{y}-\mathbf{X}\boldsymbol\beta\right)&= \arg\min_{\boldsymbol \beta}\sum_{i=1}^{n}(1/w_i)^{-1}(y_i-\mathbf{x}^{T}_i\boldsymbol\beta)^2 \\ &= \arg\min_{\boldsymbol \beta}\sum_{i=1}^{n}w_i(y_i-\mathbf{x}^{T}_i\boldsymbol\beta)^2 \end{align}$$ dan diperkirakan menggunakan dari GLS asumsi .$\boldsymbol\beta$ $$\hat{\boldsymbol\beta} = (\mathbf{X}^{T}\mathbf{V}^{-1}\mathbf{X})^{-1}\mathbf{X}^{T}\mathbf{V}^{-1}\mathbf{y}$$