Gibt es eine Möglichkeit, die Kovarianzmatrix zu verwenden, um Koeffizienten für die multiple Regression zu finden?

Für eine einfache lineare Regression kann der Regressionskoeffizient direkt aus der Varianz-Kovarianz-Matrix berechnet werden , und zwar durch wobei der Index der abhängigen Variablen und der Index der erklärenden Variablen ist.C d , $C$ d

Wenn man nur die Kovarianzmatrix hat, ist es möglich, die Koeffizienten für ein Modell mit mehreren erklärenden Variablen zu berechnen?

ETA: Für zwei erklärende Variablen scheint es, dass und analog für . Ich sehe nicht sofort, wie man dies auf drei oder mehr Variablen erweitert. β

Ja, die Kovarianzmatrix aller Variablen - erklärend und antwortend - enthält die Informationen, die zum Auffinden aller Koeffizienten erforderlich sind, vorausgesetzt, das Modell enthält einen (konstanten) Intercept-Term. (Obwohl die Kovarianzen keine Informationen über den konstanten Term liefern, können sie den Daten entnommen werden.)

Analyse

Die Daten für die erklärenden Variablen seien als dimensionale Spaltenvektoren und die Antwortvariable als Spaltenvektor , was als Realisierung einer Zufallsvariablen . Das gewöhnliche kleinste Quadrat schätzt der Koeffizienten im Modellx 1 , x 2 , ... , x p y Y $n$$x_1, x_2, \ldots, x_p$$y$$Y$$\hat\beta$

werden erhalten, indem die Spaltenvektoren zu einem Array und das System linearer Gleichungen gelöst werdenX 0 = ( 1 , 1 , … , 1 ) ' , X 1 , … , X p n × p + 1 $p+1$$X_0 = (1, 1, \ldots, 1)^\prime, X_1, \ldots, X_p$$n \times p+1$$X$

Es entspricht dem System

Die Gaußsche Eliminierung löst dieses System. Es wird fortgefahren, indem die Matrix und der -Vektor zu einem Array und zeilenreduzierend. $p+1\times p+1$p+1$\frac{1}{n}X^\prime X$$p+1$p+1×p+2$\frac{1}{n}X^\prime y$$p+1 \times p+2$$A$

Im ersten Schritt wird . Wenn festgestellt wird, dass dies ungleich Null ist, subtrahiert es geeignete Vielfache der ersten Zeile von von den verbleibenden Zeilen, um die verbleibenden Einträge in seiner ersten Spalte auf Null zu setzen. Diese Vielfachen sind und die vom Eintrag subtrahierte Zahl ist gleich . Dies ist nur die Formel für die Kovarianz von und . Darüber hinaus ist die an der Position verbleibende Zahl gleichA$\frac{1}{n}(X^\prime X)_{11} = \frac{1}{n}X_0^\prime X_0 = 1$$A$Ai+1,j+1=X ' i Xj ¯ X i ¯ X jXiXj$\frac{1}{n}X_0^\prime X_i = \overline X_i$$A_{i+1,j+1} = X_i^\prime X_j$$\overline X_i \overline X_j$$X_i$$X_j$$i+1, p+2$$\frac{1}{n}X_i^\prime y - \overline{X_i}\overline{y}$, die Kovarianz von mit .$X_i$$y$

Somit wird das System nach dem ersten Schritt der Gaußschen Eliminierung auf das Lösen reduziert

und offensichtlich - da alle Koeffizienten Kovarianzen sind - kann diese Lösung aus der Kovarianzmatrix aller Variablen gefunden werden.

(Wenn invertierbar ist, kann die Lösung geschrieben werden: . Die in der Frage angegebenen Formeln sind Spezialfälle davon, wenn und explizite Schreiben solcher Formeln wird mit wachsendem immer komplexer . Außerdem sind sie für die numerische Berechnung unterlegen, die am besten durch Lösen des Gleichungssystems und nicht durch Invertieren der Matrix .)C - 1 ( Cov ( X i , y ) ) ' p = 1 p = 2 p $C$$C^{-1}(\text{Cov}(X_i, y))^\prime$$p=1$$p=2$$p$$C$

Der konstante Term ist die Differenz zwischen dem Mittelwert von und den aus den Schätzungen vorhergesagten Mittelwerten .X $y$$X\hat{\beta}$

Beispiel

Zur Veranschaulichung werden im folgenden RCode einige Daten erstellt, deren Kovarianzen berechnet und die Koeffizientenschätzungen für die kleinsten Quadrate nur anhand dieser Informationen ermittelt. Sie werden mit den Schätzungen verglichen, die vom Schätzer der kleinsten Quadrate erhalten wurden lm.

#
# 1. Generate some data.
#
n <- 10        # Data set size
p <- 2         # Number of regressors
set.seed(17)
z <- matrix(rnorm(n*(p+1)), nrow=n, dimnames=list(NULL, paste0("x", 1:(p+1))))
y <- z[, p+1]
x <- z[, -(p+1), drop=FALSE]; 
#
# 2. Find the OLS coefficients from the covariances only.
#
a <- cov(x)
b <- cov(x,y)
beta.hat <- solve(a, b)[, 1]  # Coefficients from the covariance matrix
#
# 2a. Find the intercept from the means and coefficients.
#
y.bar <- mean(y)
x.bar <- colMeans(x)
intercept <- y.bar - x.bar %*% beta.hat  

Die Ausgabe zeigt Übereinstimmung zwischen den beiden Methoden:

(rbind(`From covariances` = c(`(Intercept)`=intercept, beta.hat),
       `From data via OLS` = coef(lm(y ~ x))))

                  (Intercept)        x1        x2
From covariances     0.946155 -0.424551 -1.006675
From data via OLS    0.946155 -0.424551 -1.006675