Die Theorie hinter dem Gewichtungsargument in R bei Verwendung von lm ()

Nach einem Jahr in der Graduiertenschule verstehe ich die "gewichteten kleinsten Quadrate" wie folgt: Sei , eine Entwurfsmatrix, \ boldsymbol \ beta \ in \ mathbb {R} ^ p ist ein Parametervektor, \ boldsymbol \ epsilon \ in \ mathbb {R} ^ n ist ein Fehlervektor, so dass \ boldsymbol \ epsilon \ sim \ mathcal {N} (\ mathbf {0}, \ sigma ^ 2 \ mathbf {V}) , wobei \ mathbf {V} = \ text {diag} (v_1, v_2, \ dots, v_n) und \ sigma ^ 2> 0 . Dann das Modell \ mathbf {y} = \ mathbf {X} \ boldsymbol \ beta + \ boldsymbol \ epsilon$\mathbf{y} \in \mathbb{R}^n$$\mathbf{X}$$n \times p$$\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$$ unter den Annahmen wird das Modell "gewichtete kleinste Quadrate" genannt. Das WLS-Problem besteht darin, $$\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}$$ Angenommen, $\mathbf{y} = \begin{bmatrix} y_1 & \dots & y_n\end{bmatrix}^{T}$ , $\boldsymbol\beta = \begin{bmatrix} \beta_1 & \dots & \beta_p\end{bmatrix}^{T}$ und $$\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$ , also $$\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{.}$$ Dies ergibt $$\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} ergibt also $$\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$ wird unter Verwendung von $$\hat{\boldsymbol\beta} = (\mathbf{X}^{T}\mathbf{V}^{-1}\mathbf{X})^{-1}\mathbf{X}^{T}\mathbf{V}^{-1}\mathbf{y}\text{.}$$ Dies ist der Umfang des Wissens, mit dem ich vertraut bin. Mir wurde nie beigebracht, wie $v_1, v_2, \dots, v_n$ sollten, obwohl es nach hier zu urteilen scheint , dass normalerweise $\text{Var}(\boldsymbol\epsilon) = \text{diag}(\sigma^2_1, \sigma^2_2, \dots, \sigma^2_n)$, was intuitiv Sinn macht. (Geben Sie stark variablen Gewichten weniger Gewicht im WLS-Problem und geben Sie Beobachtungen mit weniger Variabilität mehr Gewicht.)

Was mich besonders interessiert, ist, wie RGewichte in der lm()Funktion behandelt werden, wenn Gewichte Ganzzahlen zugewiesen werden. Von der Verwendung ?lm:

Nichtgewichte NULLkönnen verwendet werden, um anzuzeigen, dass unterschiedliche Beobachtungen unterschiedliche Abweichungen aufweisen (wobei die Werte in Gewichten umgekehrt proportional zu den Abweichungen sind); oder äquivalent, wenn die Elemente von Gewichten positive ganze Zahlen sind, ist jede Antwort der Mittelwert von Einheitsgewichtsbeobachtungen (einschließlich des Falls, dass es Beobachtungen gibt, die gleich und die Daten zusammengefasst wurden).$w_i$$y_i$$w_i$$w_i$$y_i$

Ich habe diesen Absatz mehrmals gelesen und es macht für mich keinen Sinn. Angenommen, ich habe unter Verwendung des oben entwickelten Frameworks die folgenden simulierten Werte:

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  

Wie werden diese Parameter unter Verwendung des oben entwickelten Frameworks abgeleitet? Hier ist mein Versuch, dies von Hand zu tun: Unter der Annahme von haben wir und dies in ergibt (beachten Sie, dass die Invertierbarkeit in diesem Fall nicht funktioniert, daher habe ich eine verallgemeinerte Umkehrung verwendet):$\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

Diese stimmen nicht mit den Werten aus der lm()Ausgabe überein . Was mache ich falsch?

Die Matrix sollte nicht Auch sollten Sie nicht sein .$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

Um dies genauer zu beantworten, geht die gewichtete Regression der kleinsten Quadrate unter Verwendung von weightsin von Rfolgenden Annahmen aus: Nehmen wir an, wir haben weights = c(w_1, w_2, ..., w_n). Sei , eine Entwurfsmatrix, ein Parametervektor und ist ein Fehlervektor mit dem Mittelwert und der Varianzmatrix , wobei . Dann ist Nach den gleichen Schritten der Ableitung im ursprünglichen Beitrag haben wir $\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}$$ und werden unter Verwendung von geschätzt vom GLS Annahmen .$\boldsymbol\beta$ $$\hat{\boldsymbol\beta} = (\mathbf{X}^{T}\mathbf{V}^{-1}\mathbf{X})^{-1}\mathbf{X}^{T}\mathbf{V}^{-1}\mathbf{y}$$