Мне нужно создать серию из N случайных двоичных переменных с заданной корреляционной функцией. Пусть x = {x i} - это серия двоичных переменных (с использованием значения 0 или 1, я от 1 до N). Маргинальной вероятности задается Pr (x i= 1) = p, а переменные должны быть скоррелированы следующим образом:
Corr [x i x j] = const & times; | я & minus; j | & minus; & alpha; (для я!= j)
где & alpha; - положительное число.
Если это проще, рассмотрим корреляционную функцию:
Corr [x i x j] = (| я & minus; j | +1) & minus; & alpha;
Существенная часть заключается в том, что я хочу исследовать поведение, когда корреляционная функция идет как степенной закон. (не & alpha; | я & minus; j |)
Можно ли создать такую серию, предпочтительно в Python?
Спасибо за все ваши материалы. Я нашел ответ на свой вопрос в симпатичной маленькой статье Chul Gyu Park и др., Поэтому, если кто-то столкнется с той же проблемой, посмотрите:
"Простой способ генерации коррелированных двоичных переменных" (jstor.org.stable/2684925)
для простого алгоритма. Алгоритм работает, если все элементы в корреляционной матрице положительны и для общего маргинального распределения Pr (x_i) = p_i.
J
Вы описываете случайный процесс, и мне кажется, что это сложно для меня... если вы устранили требование бинарного (0,1) и вместо этого указали ожидаемое значение и дисперсию, можно было бы описать это как генератор белого шума, питающийся через однополюсный фильтр нижних частот, который, как я думаю, даст вам характеристику & alpha; | ij |.
На самом деле это может удовлетворить бар для mathoverflow.net, в зависимости от того, как это сформулировано. Позвольте мне спросить....
update: я сделал запрос на mathoverflow.net для случая & alpha; | i-j |. Но, возможно, есть некоторые идеи, которые могут быть адаптированы к вашему делу.
Быстрый поиск по RSeek показывает, что R имеет пакеты
чтобы сделать это.
Выразите распределение x i как линейную комбинацию некоторых независимых базисных распределений f j: x i= a i1 f 1 + a i2 f 2 +.... Ограничим f j независимыми переменными, равномерно распределенными в 0..1 или в {0,1} (дискретными). Выделим теперь все, что мы знаем в матричной форме:
Let X be the vector (x1, x2, .., xn)
Let A be the matrix (a_ij) of dimension (k,n) (n rows, k columns)
Let F be the vector (f1, f2, .., fk)
Let P be the vector (p1, p2, .., pn)
Let R be the matrix (E[x_i,x_j]) for i,j=1..n
Definition of the X distribution: X = A * F
Constraint on the mean of individual X variables: P = A * (1 ..k times.. 1)
Correlation constraint: AT*A = 3R or 2R in the discrete case (because E[x_i x_j] =
E[(a_i1*f_1 + a_i2*f_2 + ...)*(a_j1*f_1 + a_j2*f_2 + ...)] =
E[sum over p,q: a_ip*f_p*a_jq*f_q] = (since for p/=q holds E[f_p*f_q]=0)
E[sum over p: a_ip*a_jp*f_p^2] =
sum over p: a_ip*a_jp*E[f_p^2] = (since E[f_p^2] = 1/3 or 1/2 for the discrete case)
sum over p: 1/3 or 1/2*a_ip*a_jp
And the vector consisting of those sums over p: a_ip*a_jp is precisely AT*A.
Теперь вам нужно решить два уравнения:
AT*A = 3R (or 2R in the discrete case)
A*(1...1) = P
Я прошу математиков исправить меня, потому что я вполне мог бы смешать ATA с AAT или сделать что-то еще более неправильное.
Чтобы сгенерировать значение x i как линейную смесь базисных распределений, используйте двухэтапный процесс: 1) используйте единую случайную величину, чтобы выбрать одно из базисных распределений, взвешенное с соответствующая вероятность, 2) генерирует результат с использованием выбранного базисного распределения.
Решение грубой силы должно выражать ограничения задачи как линейной программы с переменными 2^N pr(w), где w пробегает все двоичные строки длины N. Сначала ограничение, что pr - распределение вероятности:
for all w: 0 <= pr(w) <= 1
sum_w pr(w) = 1
Во-вторых, ограничение того, что ожидание каждой переменной будет p:
for all i: sum_{w such that w[i] = 1} pr(w) = p
В-третьих, ограничения ковариации:
for all i < j: sum_{w such that w[i] = w[j] = 1} pr(w) = const * |j - i|^alpha - p^2
Это очень медленно, но беглый поиск литературы не стал лучше. Если вы решите его реализовать, вот некоторые решатели LP с привязками Python: http://wiki.python.org/moin/NumericAndScientific/Libraries
Здесь есть интуитивный/экспериментальный подход, который, кажется, работает.
Если b является двоичным r.v., m - среднее значение двоичного r.v., c - это необходимая корреляция, rand() генерирует U (0,1) r.v. и d является коррелированным двоичным r.v. вы хотите:
d = if (rand() < c, b, if (rand() < m, 0, 1))
То есть, если равномерное r.v. меньше желаемой корреляции, d = b. В противном случае d = другое случайное двоичное число.
Я запускал этот 1000 раз для столбца 2000 двоичных файлов r.v.s. с m =.5 и c =.4 и c =.5 Среднее значение корреляции было точно таким же, как указано, распределение оказалось нормальным. Для корреляции 0,4 std отклонение корреляции составляло 0,02.
Извините - я не могу доказать, что это работает все время, но вы должны признать, что это легко.