Я не могу понять, как создать многошаговый код с двумерным нормальным распределением

Question

Рассмотрим следующий метод для генерации случайной выборки из двумерного нормального распределения. Сначала сгенерируйте случайное число x, y независимо от стандартной нормали и сформируйте пару (x, y) ^ '; во-вторых, вычислите факторизацию Холецкого, где ρ - некоторая постоянная, принимает значение в (−1, 1), а R - матрица верхнего треугольника; в-третьих, вычислите новую пару. Получившаяся новая пара имеет двумерное нормальное распределение с коэффициентом корреляции ρ. затем сгенерировать 100 пар случайных векторов с ρ = −0,8, 0,5, 0,8 соответственно; затем) построить график рассеяния y против x для каждого из трех случаев. для каждого из трех случаев вычислите среднее значение выборки, выборочную дисперсию x и y соответственно и выборочную корреляцию между x и y

Answer 1

Сначала несколько критических комментариев:

1. Это хорошее упражнение, и я настоятельно рекомендую вам потратить некоторое время на изучение реализаций кодирования. Вы можете найти много соответствующих материалов в Интернете.
2. Что касается домашних заданий, мы обычно ожидаем, что вы будете придерживаться некоторых ключевых правил ; Короче говоря, мы ожидаем, что вы продемонстрируете подлинную попытку решить проблему самостоятельно. Вы должны поделиться своей попыткой кода и четко указать, где вы застряли. Никто здесь не любит вопросы "gimmeh teh codez" .

С этими вещами в стороне, вот некоторые «идеи» кода, которые должны помочь вам с правильным решением вашей домашней работы. Это не полное решение , потому что, как объяснено выше и в комментариях, это не так, как работает SO.

---

* * Фон тысяча двадцать-один * * тысяча двадцать-дв

В разложении Холецкого мы можем разложить положительно определенную матрицу дисперсии-ковариации Sigma как Sigma = R R^T для некоторой нижней треугольной матрицы R. Вектор X = mu + R Z тогда имеет многомерное нормальное распределение, где

• Z - это вектор (Z_1, Z_2, ..., Z_n) из n независимых стандартных нормальных переменных,
• R - нижняя треугольная матрица r x n, связанная с дисперсионно-ковариационной матрицей посредством разложения Холецкого, и
• mu - это вектор (mu_1, mu_2, ..., mu_r) средних.

Для стандартного нормального случая bivariate дисперсионно-ковариационная матрица представляет собой просто матрицу 2 x 2 с диагональю единицы и rho на недиагональных элементах.

R код

Мы можем определить функцию bvsigma, которая возвращает двумерную дисперсионно-ковариационную матрицу для данного коэффициента корреляции rho и дисперсий sigma_1 и sigma_2.

bvsigma <- function(rho = 0.8, sigma_1 = 1, sigma_2 = 1)
    matrix(c(sigma_1^2, rho * sigma_1 * sigma_2, rho * sigma_1 * sigma_2, sigma_2^2), ncol = 2)

Для стандартного двумерного нормального распределения с rho = 0.8 мы имеем

sigma <- bvsigma(rho = 0.8)
sigma
#     [,1] [,2]
#[1,]  1.0  0.8
#[2,]  0.8  1.0

Тогда нижняя треугольная матрица будет

R <- t(chol(sigma))
R
#     [,1] [,2]
#[1,]  1.0  0.0
#[2,]  0.8  0.6

Мы подтверждаем, что действительно R R^T восстанавливает матрицу дисперсии-ковариации

R %*% t(R)
#     [,1] [,2]
#[1,]  1.0  0.8
#[2,]  0.8  1.0

Теперь мы можем собрать все части вместе и определить функцию rbvnorm, которая возвращает n выборок из двумерного стандартного нормаля с матрицей дисперсии-ковариации sigma

rbvnorm <- function(n, sigma) {
    R <- t(chol(sigma))
    t(do.call(cbind, replicate(n, {
        Z <- rnorm(2)
        R %*% Z
    }, simplify = FALSE)))}

Давайте нарисуем n = 1000 выборок из двумерного стандартного нормального распределения с rho = 0.8 и нанесем данные

set.seed(2018)
mat <- rbvnorm(1000, bvsigma(rho = 0.8))

# Plot
library(ggplot2)
ggplot(data = as.data.frame(mat), aes(V1, V2)) +
    geom_point()