La teoría detrás del argumento de los pesos en R cuando se usa lm ()

Después de un año en la escuela de posgrado, mi comprensión de los "mínimos cuadrados ponderados" es la siguiente: let , sea ​​una matriz de diseño , \ boldsymbol \ beta \ in \ mathbb {R} ^ p sea ​​un vector de parámetros, \ boldsymbol \ epsilon \ in \ mathbb {R} ^ n sea ​​un vector de error tal que \ boldsymbol \ epsilon \ sim \ mathcal {N} (\ mathbf {0}, \ sigma ^ 2 \ mathbf {V}) , donde \ mathbf {V} = \ text {diag} (v_1, v_2, \ dots, v_n) y \ sigma ^ 2> 0 . Entonces el modelo \ 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$ ϵ ∼ N ( 0 , σ 2 V$\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$$ bajo los supuestos se llama el modelo de "mínimos cuadrados ponderados". El problema de WLS termina siendo encontrar $$\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}$$ Supongamos que $\mathbf{y} = \begin{bmatrix} y_1 & \dots & y_n\end{bmatrix}^{T}$ , $\boldsymbol\beta = \begin{bmatrix} \beta_1 & \dots & \beta_p\end{bmatrix}^{T}$ y $$\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$ , entonces $$\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{.}$$ Esto da $$\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} dando así $$\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$ se estima usando $$\hat{\boldsymbol\beta} = (\mathbf{X}^{T}\mathbf{V}^{-1}\mathbf{X})^{-1}\mathbf{X}^{T}\mathbf{V}^{-1}\mathbf{y}\text{.}$$ Esta es la extensión del conocimiento con el que estoy familiarizado. Nunca me enseñaron cómo deberían elegirse $v_1, v_2, \dots, v_n$ , aunque parece que, a juzgar por esto , generalmente $\text{Var}(\boldsymbol\epsilon) = \text{diag}(\sigma^2_1, \sigma^2_2, \dots, \sigma^2_n)$, lo que tiene sentido intuitivo. (Proporcione pesos muy variables menos peso en el problema WLS, y brinde observaciones con menos variabilidad más peso).

Lo que me interesa especialmente es cómo Rmaneja los pesos en la lm()función cuando los pesos se asignan como enteros. De usar ?lm:

Las no NULLponderaciones se pueden usar para indicar que las diferentes observaciones tienen diferentes variaciones (con los valores en pesos inversamente proporcionales a las variaciones); o de manera equivalente, cuando los elementos de los pesos son enteros positivos , que cada respuesta es la media de las observaciones de peso unitario (incluido el caso de que hay observaciones iguales a y los datos se han resumido).$w_i$$y_i$$w_i$$w_i$$y_i$

He releído este párrafo varias veces, y no tiene sentido para mí. Usando el marco que desarrollé anteriormente, supongamos que tengo los siguientes valores simulados:

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  

Usando el marco que he desarrollado anteriormente, ¿cómo se derivan estos parámetros? Aquí está mi intento de hacer esto a mano: suponiendo , tenemos y hacer esto en give (tenga en cuenta que la invertibilidad no funciona en este caso, por lo que utilicé un inverso generalizado):$\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

Estos no coinciden con los valores de la lm()salida. ¿Qué estoy haciendo mal?

La matriz debe ser no Además, tu deberías ser , no .$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

Para responder esto de manera más concisa, la regresión ponderada de mínimos cuadrados usando weightsin Rhace los siguientes supuestos: supongamos que tenemos weights = c(w_1, w_2, ..., w_n). Deje que , sea ​​una matriz de diseño , sea ​​un vector de parámetros y puede ser un vector de error con media y varianza matriz , donde . Luego, Siguiendo los mismos pasos de la derivación en la publicación original, tenemos $\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}$$ y se estima usando del GLS supuestos .$\boldsymbol\beta$ $$\hat{\boldsymbol\beta} = (\mathbf{X}^{T}\mathbf{V}^{-1}\mathbf{X})^{-1}\mathbf{X}^{T}\mathbf{V}^{-1}\mathbf{y}$$