scipy - генерировать случайные величины с корреляциями

Question

Я работаю над реализацией базового симулятора Монте-Карло в Python для некоторого моделирования рисков управления проектами, который я пытаюсь сделать (в основном Crystal Ball/@Risk, но в Python).

У меня есть набор из n случайных переменных (все scipy.stats экземпляров). Я знаю, что могу использовать rv.rvs(size=k) для генерации kнезависимых наблюдений от каждой из этих n переменных.

Я хотел бы ввести корреляции между переменными, указав положительную полуопределенную корреляционную матрицу n x n.

Есть ли чистый способ сделать это в scipy?

Что я Пробовал

This answer и this answer, кажется, указывают, что «копулы» будет ответ, но я не вижу каких-либо ссылок на SciPy к ним.

This link, похоже, реализует то, что я ищу, но я не уверен, что эта функция реализована уже в scipy. Мне также хотелось бы, чтобы он работал для ненормальных переменных.

Кажется, что стандартный метод Iman, Conover paper.

Answer 1

Я просто хочу корреляцию через Gaussian Copula (*), тогда ее можно вычислить в несколько шагов с помощью numpy и scipy.

• создавать многомерные случайные величины с требуемой ковариации, numpy.random.multivariate_normal и создание (Ноббс по k_variables) массива

• применять scipy.stats.norm.cdf для преобразования нормально равномерных случайных величин, для каждого столбца/переменной, чтобы получить равномерное маргинальный распределения

• dist.ppf применяется для преобразования равномерной маржи до требуемого распределения, где dist может быть одним из распределений в scipy.stats

(*) Gaussian копула только один выбор, и это не самое лучшее, когда мы заинтересованы в хвостовом поведении, но это проще всего работать с , например http://archive.wired.com/techbiz/it/magazine/17-03/wp_quant?currentPage=all

две ссылок

https://stats.stackexchange.com/questions/37424/how-to-simulate-from-a-gaussian-copula

http://www.mathworks.com/products/demos/statistics/copulademo.html

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

Answer 2

Похоже, что метод выборки на основе отбраковки, такой как алгоритм Metropolis-Hastings, - это то, что вы хотите. Scipy может реализовать такие методы с помощью функции scipy.optimize.basinhopping.

Методы выборки на основе отклонений позволяют вам отбирать образцы из любого заданного распределения вероятностей. Идея состоит в том, что вы производите произвольные выборки из другого «предложения» pdf, который легко пробовать (например, для унифицированных или гауссовых дистрибутивов), а затем использовать случайный тест, чтобы решить, должен ли этот образец из распределения предложений быть «принят» как представляющий образец желаемого распределения.

Оставшихся уловок тогда будет:

1. Выяснить формы совместной N-мерной функции плотности вероятности, которая имеет маргинал формы вы хотите по каждому измерению, но с корреляционной матрицей, что вы хотеть. Это легко сделать для распределения Гаусса, где требуемая корреляционная матрица и средний вектор - это все, что вам нужно для определения распределения. Если у ваших маргиналов есть простое выражение, вы можете найти этот pdf-файл с помощью простой, но утомительной алгебры. This бумага цитирует несколько других, которые делают то, о чем вы говорите, и я уверен, что их гораздо больше.

2. Сформулируйте функцию для basinhopping, чтобы свести к минимуму так, чтобы принятое количество «минимумов» составляло образцы этого pdf, которые вы определили.

Учитывая результаты (1), (2), должно быть просто.