Глава 3: The Multivariate Normal Distribution

[4.2]

Да се покаже дека \(\boldsymbol{z}=(\boldsymbol{T'})^{-1}(\boldsymbol{y}-\boldsymbol{\mu}) \implies E(\boldsymbol{z})=\boldsymbol{0},cov(\boldsymbol{z})=\boldsymbol{I}\)

\[ E(\boldsymbol{z})=E((\boldsymbol{T'})^{-1}(\boldsymbol{y}-\boldsymbol{\mu}))=(\boldsymbol{T'})^{-1}E(\boldsymbol{y}-\boldsymbol{\mu})=(\boldsymbol{T'})^{-1}(E(\boldsymbol{y})-\boldsymbol{\mu})=(\boldsymbol{T'})^{-1}(\boldsymbol{\mu}-\boldsymbol{\mu})=\boldsymbol{0}\\ cov(\boldsymbol{z})=cov((\boldsymbol{T'})^{-1}(\boldsymbol{y}-\boldsymbol{\mu}))= (\boldsymbol{T}')^{-1}cov(\boldsymbol{y})[(\boldsymbol{T}')^{-1}]'= (\boldsymbol{T}')^{-1}\boldsymbol{T}'\boldsymbol{T}\boldsymbol{T}^{-1}=\boldsymbol{I} \]

[4.5]

Да се покаже дека \[\frac{1}{2}\sum_{i=1}^{n}(\boldsymbol{y}_i-\bar{\boldsymbol{y}}+\bar{\boldsymbol{y}}-\boldsymbol{\mu})'\boldsymbol{\Sigma}^{-1}(\boldsymbol{y}_i-\bar{\boldsymbol{y}}+\bar{\boldsymbol{y}}-\boldsymbol{\mu})=\frac{1}{2}\sum_{i=1}^{n}(\boldsymbol{y}_i-\bar{\boldsymbol{y}})'\boldsymbol{\Sigma}^{-1}(\boldsymbol{y}_i-\bar{\boldsymbol{y}})+\frac{n}{2}(\bar{\boldsymbol{y}}-\boldsymbol{\mu})'\boldsymbol{\Sigma}^{-1}(\bar{\boldsymbol{y}}-\boldsymbol{\mu})\]

\[ \frac{1}{2}\sum_{i=1}^{n}(\boldsymbol{y}_i-\bar{\boldsymbol{y}}+\bar{\boldsymbol{y}}-\boldsymbol{\mu})'\boldsymbol{\Sigma}^{-1}(\boldsymbol{y}_i-\bar{\boldsymbol{y}}+\bar{\boldsymbol{y}}-\boldsymbol{\mu})=\\ \frac{1}{2}\sum_{i=1}^{n}(\boldsymbol{y}_i-\bar{\boldsymbol{y}})'\boldsymbol{\Sigma}^{-1}(\boldsymbol{y}_i-\bar{\boldsymbol{y}})+ \frac{1}{2}\sum_{i=1}^{n}(\boldsymbol{y}_i-\bar{\boldsymbol{y}})'\boldsymbol{\Sigma}^{-1}(\bar{\boldsymbol{y}}-\boldsymbol{\mu})+ \frac{1}{2}\sum_{i=1}^{n}(\bar{\boldsymbol{y}}-\boldsymbol{\mu})'\boldsymbol{\Sigma}^{-1}(\boldsymbol{y}_i-\bar{\boldsymbol{y}})+ \frac{1}{2}\sum_{i=1}^{n}(\bar{\boldsymbol{y}}-\boldsymbol{\mu})'\boldsymbol{\Sigma}^{-1}(\bar{\boldsymbol{y}}-\boldsymbol{\mu})=\\ \frac{1}{2}\sum_{i=1}^{n}(\boldsymbol{y}_i-\bar{\boldsymbol{y}})'\boldsymbol{\Sigma}^{-1}(\boldsymbol{y}_i-\bar{\boldsymbol{y}})+ \frac{1}{2}(n\bar{\boldsymbol{y}}-n\bar{\boldsymbol{y}})'\boldsymbol{\Sigma}^{-1}(\bar{\boldsymbol{y}}-\boldsymbol{\mu})+ \frac{1}{2}(\bar{\boldsymbol{y}}-\boldsymbol{\mu})'\boldsymbol{\Sigma}^{-1} (n\bar{\boldsymbol{y}}-n\bar{\boldsymbol{y}})+ \frac{n}{2}(\bar{\boldsymbol{y}}-\boldsymbol{\mu})'\boldsymbol{\Sigma}^{-1}(\bar{\boldsymbol{y}}-\boldsymbol{\mu})= \\ \frac{1}{2}\sum_{i=1}^{n}(\boldsymbol{y}_i-\bar{\boldsymbol{y}})'\boldsymbol{\Sigma}^{-1}(\boldsymbol{y}_i-\bar{\boldsymbol{y}})+\frac{n}{2}(\bar{\boldsymbol{y}}-\boldsymbol{\mu})'\boldsymbol{\Sigma}^{-1}(\bar{\boldsymbol{y}}-\boldsymbol{\mu}) \]

[4.9]

Да се покаже дека \(F_n=\frac{n-p-1}{p}(\frac{1}{w}-1)\)

\[ w=1-\frac{nD_n^2}{(n-1)^2} \implies F_n=\frac{n-p-1}{p}(\frac{1}{1-\frac{nD_n^2}{(n-1)^2}}-1)=\frac{n-p-1}{p}(\frac{1}{w}-1) \]

[4.11]

\[ \boldsymbol{y} \sim N_3(\boldsymbol{\mu}, \boldsymbol{\Sigma}); \boldsymbol{\mu}=\begin{pmatrix}3\\1\\4\end{pmatrix}; \boldsymbol{\Sigma}=\begin{pmatrix}6 & 1 & -2\\1 & 13 & 4\\-2 & 4 & 4\end{pmatrix} \]

(a)

Да се најде \(\boldsymbol{z}=(\boldsymbol{T}')^{-1}(\boldsymbol{y}-\boldsymbol{\mu}) \sim N_3(\boldsymbol{0},\boldsymbol{I})\)

Чолески декомпозиција:

cholesky <- function(a) {
  d = dim(a)[1]
  t = matrix(rep(0, d * d), nrow=d, ncol=d)
  t[1, 1] = sqrt(a[1, 1])
  for (j in 2:d) {
    t[1, j] = a[1, j] / t[1, 1]
  }
  for (i in 2:d) {
    for (j in 2:d) {
      s = 0
      if (i == j) {
        for (k in 1:(i-1)) {
          s = s + t[k, i] ^ 2
        }
        t[i, i] = sqrt(a[i, i] - s)
      } else if (i < j) {
        for (k in 1:(i-1)) {
          s = s + t[k, i] * t[k, j]
        }
        t[i, j] = (a[i, j] - s) / t[i, i]
      }
    }
  }
  return(t)
}

Sigma = matrix(
  c(6, 1, -2,
    1, 13, 4,
    -2, 4, 4
  ),
  nrow=3,
  ncol=3
)
T = cholesky(Sigma)
T

        [,1]      [,2]       [,3]
[1,] 2.44949 0.4082483 -0.8164966
[2,] 0.00000 3.5823642  1.2096295
[3,] 0.00000 0.0000000  1.3675269

Векторот \(\boldsymbol{z}\):

\[ \boldsymbol{z} = \begin{pmatrix} 2.45 & 00 & 0\\ 0.41 & 3.58 & 0\\ -0.82 & 1.21 & 1.37\end{pmatrix}^{-1}(\boldsymbol{y}-\boldsymbol{\mu}) \]

(b)

Да се најде \(\boldsymbol{z}=(\boldsymbol{\Sigma}^{\frac{1}{2}})^{-1}(\boldsymbol{y}-\boldsymbol{\mu}) \sim N_3(\boldsymbol{0},\boldsymbol{I})\)

Прво правиме спектрална декомпозиција на \(\boldsymbol{\Sigma}\). Матрицата е симетрична, затоа има 3 сопствени вредности. Карактеристичната равенка изгледа вака:

\[ |\boldsymbol{\Sigma}-\lambda\boldsymbol{I}|=\begin{vmatrix}6-\lambda & 1 & -2\\1 & 13-\lambda & 4\\-2 & 4 & 4-\lambda\end{vmatrix}=\\ (6-\lambda)\begin{vmatrix}13-\lambda & 4\\4 & 4-\lambda\end{vmatrix} -\begin{vmatrix}1 & 4\\-2 & 4-\lambda\end{vmatrix} -2\begin{vmatrix}1 & 13-\lambda\\-2 & 4\end{vmatrix}=\\ (6-\lambda)(36-17\lambda+\lambda^2)-(12-\lambda)-2(30-2\lambda)=\\ 144-133\lambda+23\lambda^2-\lambda^3 \]

Корените на равенката се:

lambdas = polyroot(z=c(144, -133, 23, -1))
print(lambdas)

[1]  1.401827+0i  7.071196-0i 14.526977+0i

За да ги најдеме сопствените вектори, потребно ни е да го решиме системот равенки:

\[ \begin{pmatrix}6-\lambda & 1 & -2\\1 & 13-\lambda & 4\\-2 & 4 & 4-\lambda\end{pmatrix} \begin{pmatrix}x_1\\x_2\\x_3\end{pmatrix}= \begin{pmatrix}0\\0\\0\end{pmatrix} \]

Во овој систем знаеме дека едната равенка е вишок. Затоа, претпоставувајќи дека \(x_3=1\), додавајќи ја последната равенка на првите две го добиваме системот:

\[ \begin{pmatrix} (6-\lambda)-2 & (1)+4\\ (1)-2 & (13-\lambda)+4 \end{pmatrix} \begin{pmatrix}x_1\\x_2\end{pmatrix}= \begin{pmatrix} -[(-2) + (4 - \lambda)]\\ -[(4)+(4-\lambda)] \end{pmatrix} \]

Овој систем го решаваме со кодот од домашната за Глава 2:

solve_equations <- function(a, b) {
  a <- matrix(data=a, nrow=2, ncol=2, byrow=TRUE)
  b <- matrix(data=b, nrow=2, ncol=1, byrow=FALSE)
  sol <- solve(a, b)
  nu <- c(sol[1], sol[2], 1)
  nu_length <- sqrt(sum(nu ^ 2))
  return(nu / nu_length)
}

solve_nu <- function(l) {
  a <- c(
    6 - l - 2, 1 + 4, 
    1 - 2, 13 - l + 4)
  b <- c(
    2 - (4 - l), 
    -4 - (4 - l)
  )
  return(solve_equations(a, b))
}

Сопствениoт вектор за \(\lambda_1\) со модул 1 e:

nu_1 = solve_nu(lambdas[1])
print(nu_1)

[1]  0.4358269+0i -0.3267970+0i  0.8386052-0i

Сопствениoт вектор за \(\lambda_2\) со модул 1 e:

nu_2 = solve_nu(lambdas[2])
print(nu_2)

[1] -0.8996077-0i -0.1296166-0i  0.4170197-0i

Сопствениoт вектор за \(\lambda_3\) со модул 1 e:

nu_3 = solve_nu(lambdas[3])
print(nu_3)

[1] 0.02758366+0i 0.93616414+0i 0.35047945-0i

Конечно ги добиваме \(C\) и \(D\) (оцде се користи %*% за множење матрици со цел да се провери дека \(CDC'\) е исто со почетната матрица):

C = matrix(
  data=as.numeric(c(nu_1, nu_2, nu_3)), 
  nrow=3, 
  ncol=3, 
  byrow=FALSE
)

imaginary parts discarded in coercion

D = diag(x=as.numeric(lambdas), nrow=3, ncol=3)

imaginary parts discarded in coercion

C %*% D %*% t(C)

     [,1] [,2] [,3]
[1,]    6    1   -2
[2,]    1   13    4
[3,]   -2    4    4

Сега наоѓаме \(\boldsymbol{\Sigma}^{\frac{1}{2}}=\boldsymbol{C}\boldsymbol{D}^{\frac{1}{2}}\boldsymbol{C}'\). Истовремено и проверуваме дали е точно најдено со тоа што го креваме на квадрат и наоѓаме и инверзна матрица..

D_sr = diag(x=as.numeric(sqrt(lambdas)), nrow=3, ncol=3)

imaginary parts discarded in coercion

Sigma_sr = C %*% D_sr %*% t(C)
Sigma_sr

           [,1]      [,2]       [,3]
[1,]  2.3798448 0.2398604 -0.5280208
[2,]  0.2398604 3.5114681  0.7823418
[3,] -0.5280208 0.7823418  1.7632740

Sigma_sr %*% Sigma_sr

     [,1] [,2] [,3]
[1,]    6    1   -2
[2,]    1   13    4
[3,]   -2    4    4

solve(Sigma_sr)

            [,1]        [,2]       [,3]
[1,]  0.46496849 -0.06966936  0.1701484
[2,] -0.06966936  0.32645938 -0.1657086
[3,]  0.17014841 -0.16570861  0.6916013

Решението е:

\[ \boldsymbol{z}=(\boldsymbol{\Sigma}^{\frac{1}{2}})^{-1}(\boldsymbol{y}-\boldsymbol{\mu})=\begin{pmatrix} 0.46 & -0.07 & 0.17\\ -0.07 & 0.33 & -0.17\\ 0.17 & -0.17 & 0.69 \end{pmatrix}(\boldsymbol{y}-\boldsymbol{\mu}) \sim N_3(\boldsymbol{0},\boldsymbol{I}) \]

(c)

Дистрибуцијата на \((\boldsymbol{y}-\boldsymbol{\mu})'\boldsymbol{\Sigma}^{-1}(\boldsymbol{y}-\boldsymbol{\mu})\) е хи-квадратна со 3 степени на слобода, т.е. \((\boldsymbol{y}-\boldsymbol{\mu})'\boldsymbol{\Sigma}^{-1}(\boldsymbol{y}-\boldsymbol{\mu}) \sim \chi^2_3\).

[4.13]

\[ \boldsymbol{y} \sim N_4(\boldsymbol{\mu}, \boldsymbol{\Sigma}); \boldsymbol{\mu}=\begin{pmatrix}-2\\3\\-1\\5\end{pmatrix}; \boldsymbol{\Sigma}=\begin{pmatrix} 11 & −8 & 3 & 9\\ −8 & 9 & −3 & 6\\ 3 & −3 & 2 & 3\\ 9 & 6 & 3 & 9 \end{pmatrix} \]

(a)

Да се најде \(\boldsymbol{z}=(\boldsymbol{T}')^{-1}(\boldsymbol{y}-\boldsymbol{\mu}) \sim N_4(\boldsymbol{0},\boldsymbol{I})\)

Чолески декомпозиција со кодот од претходната задача:

Sigma = matrix(
  c(11, -8, 3, 9, 
    -8, 9, -3, 6, 
    3, -3, 2, 3, 
    9, 6, 3, 9),
  nrow=4,
  ncol=4
)
T = cholesky(Sigma)

NaNs produced

T

         [,1]      [,2]       [,3]     [,4]
[1,] 3.316625 -2.412091  0.9045340 2.713602
[2,] 0.000000  1.783765 -0.4586825 7.033131
[3,] 0.000000  0.000000  0.9856108 3.826489
[4,] 0.000000  0.000000  0.0000000      NaN

Една од вредностите е NaN бидејќи \(\boldsymbol{\Sigma}\) не е позитивно дефинитна (веројатно има некоја грешка во задачата).

(b)

Да се најде \(\boldsymbol{z}=(\boldsymbol{\Sigma}^{\frac{1}{2}})^{-1}(\boldsymbol{y}-\boldsymbol{\mu}) \sim N_4(\boldsymbol{0},\boldsymbol{I})\)

Не може да се најде бидејќи \(\boldsymbol{\Sigma}\) не е позитивно дефинитна (веројатно има некоја грешка во задачата).

(c)

Дистрибуцијата на \((\boldsymbol{y}-\boldsymbol{\mu})'\boldsymbol{\Sigma}^{-1}(\boldsymbol{y}-\boldsymbol{\mu})\) е хи-квадратна со 4 степени на слобода, т.е. \((\boldsymbol{y}-\boldsymbol{\mu})'\boldsymbol{\Sigma}^{-1}(\boldsymbol{y}-\boldsymbol{\mu}) \sim \chi^2_4\).