Использование R для определения корреляций для испытаний LatinHypercube / Monte Carlo

Question

В настоящее время я использую Python и RPY для использования функциональности внутри R.

Как использовать библиотеку R для генерации образцов Монте-Карло, которые учитывают корреляцию между двумя переменными .. например если переменные A и B имеют корреляцию 85% (0,85), мне нужно сгенерировать все выборки Монте-Карло, соблюдая эту корреляцию между A & B.

Буду признателен, если кто-нибудь сможет поделиться идеями / фрагментами

Спасибо

Answer 1

Метод ранговой корреляции Имана и Коновера представляется широко используемым и общим подходом к получению коррелированных образцов Монте-Карло для компьютерных экспериментов, анализа чувствительности и т. Д. К сожалению, я только что столкнулся с этим и у вас нет доступа к PDF, поэтому вы не знаете, как авторы на самом деле реализуют свой метод, но вы могли бы продолжить.

Их метод является более общим, поскольку каждая переменная может происходить из другого распределения, в отличие от многовариантной нормали ответа @ Dirk.

Обновление: Я нашел реализацию R вышеуказанного подхода в пакете mc2d, в частности, вам нужна функция cornode().

Вот пример, взятый из ?cornode

> require(mc2d)
> x1 <- rnorm(1000)
> x2 <- rnorm(1000)
> x3 <- rnorm(1000)
> mat <- cbind(x1, x2, x3)
> ## Target
> (corr <- matrix(c(1, 0.5, 0.2, 0.5, 1, 0.2, 0.2, 0.2, 1), ncol=3))
     [,1] [,2] [,3]
[1,]  1.0  0.5  0.2
[2,]  0.5  1.0  0.2
[3,]  0.2  0.2  1.0
> ## Before
> cor(mat, method="spearman")
            x1         x2          x3
x1  1.00000000 0.01218894 -0.02203357
x2  0.01218894 1.00000000  0.02298695
x3 -0.02203357 0.02298695  1.00000000
> matc <- cornode(mat, target=corr, result=TRUE)
Spearman Rank Correlation Post Function
          x1        x2        x3
x1 1.0000000 0.4515535 0.1739153
x2 0.4515535 1.0000000 0.1646381
x3 0.1739153 0.1646381 1.0000000

Ранговые корреляции в matc теперь очень близки к целевым корреляциям corr.

Идея заключается в том, что вы рисуете выборки отдельно от распределения для каждой переменной, а затем используете подход Iman & Connover, чтобы сделать выборки (как можно ближе) к целевым корреляциям, насколько это возможно.

Answer 2

Это часто задаваемые вопросы. Вот один ответ, используя рекомендуемый пакет:

R> library(MASS)
R> example(mvrnorm)

mvrnrmR> Sigma <- matrix(c(10,3,3,2),2,2)

mvrnrmR> Sigma
     [,1] [,2]
[1,]   10    3
[2,]    3    2

mvrnrmR> var(mvrnorm(n=1000, rep(0, 2), Sigma))
        [,1]    [,2]
[1,] 8.82287 2.63987
[2,] 2.63987 1.93637

mvrnrmR> var(mvrnorm(n=1000, rep(0, 2), Sigma, empirical = TRUE))
     [,1] [,2]
[1,]   10    3
[2,]    3    2
R> 

Переключение между корреляцией и ковариацией является простым (подсказка: внешнее произведение вектора стандартных отклонений).

Answer 3

Этот вопрос не был помечен как python, но, исходя из вашего комментария, похоже, что вы, возможно, ищете и решение Python. Самая базовая Python-реализация Iman Convover, которую я могу придумать, выглядит в Python следующим образом (фактически numpy):

def makeCorrelated( y,  corMatrix ):
    c = multivariate_normal(zeros(size( y, 0 ) ) , corMatrix, size( y, 1 ) )
    key = argsort( argsort(c, axis=0), axis=0 ).T
    out = map(take, map(sort, y), key)
    out = array(out)
    return out

где y - массив выборок из маргинальных распределений, а corMatrix - положительная полуопределенная симметричная корреляционная матрица. Учитывая, что эта функция использует multivariate_normal () для матрицы c, вы можете сказать, что она использует подразумеваемую гауссову копулу. Чтобы использовать разные структуры связок, вам нужно использовать разные драйверы для матрицы c.