Matlab mvnrnd в gsl

Привет, я не уверен, что мой алгоритм правильный. Я пытаюсь воспроизвести mvnrnd. функция Matlab, но в gsl. В некоторых журнальных статьях я нашел алгоритм, который создает вектор многомерных нормальных чисел, но мне нужна матрица многомерных нормальных случайных чисел.

скажем, распределение Z ~ (мю, сигма);

предполагая, что сигма - матрица, уже положительно определенная.

алгоритм, который я нашел в Интернете, говорит

1. cholskey(sigma) = A
2. generate uniform gaussian vector r
3. matrix vector triangular product with gsl_blas_dtrmv A * r
4. add mu to Ar and that will be a vector of multivariate normal random numbers

мой метод ниже

Правильно ли приведенные ниже изменения для произведения матрицы случайных величин

    1. cholskey(sigma) = A
    2. generate uniform gaussian matrix R
    3. matrix matrix scalar product AR
    4. add mu to AR and that will be a matrix of multivariate normal random numbers

person pyCthon    schedule 01.05.2012    source источник
comment
Привет, я пытаюсь сделать подобные вещи здесь. Могу ли я узнать, как вы сгенерировали равномерный вектор Гаусса в gsl. Я видел, что упоминаются методы box-muller. Но будет ли иметь значение, если я просто сгенерирую гауссово случайное число одно за другим для вектора? Будет ли ковариация из-за псевдослучайного числа.   -  person Coconut    schedule 16.05.2019


Ответы (1)


Да это верно. См., например. эта статья Википедии о многомерных нормальных ГСЧ, в которой есть этот раздел:

Получение значений из распределения

Широко используемый метод рисования случайного вектора x из N-мерного многомерного нормального распределения со средним вектором μ и ковариационной матрицей Σ работает следующим образом:

  1. Найдите любую вещественную матрицу A такую, что A AT = Σ. Когда Σ является положительно определенным, обычно используется разложение Холецкого. [...]

  2. Пусть z = (z1,...,zN)T — вектор, компонентами которого являются N независимых стандартных нормальных переменных (которые можно сгенерировать, например, с помощью преобразования Бокса–Мюллера).

  3. Пусть x будет µ + Az. Это имеет желаемое распределение из-за свойства аффинного преобразования.

который описывает тот же алгоритм.

Также существует несколько реализаций для R, например mvrnorm в MASS, который поставляется с каждой установкой R.

person Dirk Eddelbuettel    schedule 02.05.2012
comment
Благодарность! не думал, что найду такого знаменитого, как ты, чтобы он ответил на мой вопрос - person pyCthon; 02.05.2012