2015-01-01 14 views
9

Sto lavorando per implementare un simulatore Monte Carlo di base in Python per alcuni modelli di gestione del rischio di gestione che sto cercando di fare (fondamentalmente Crystal Ball/@Risk, ma in Python).scipy - genera variabili casuali con correlazioni

Ho un set di variabili casuali n (tutte le istanze scipy.stats). So che posso utilizzare rv.rvs(size=k) per generare kindipendente da osservazioni da ognuna di queste variabili n.

Vorrei introdurre le correlazioni tra le variabili specificando una matrice di correlazione semifissa positiva n x n.

C'è un modo pulito per farlo in scipy?

quello che ho provato

This answer e this answer sembrano indicare che "copule" sarebbe una risposta, ma non vedo alcun riferimento a SciPy a loro.

This link sembra implementare quello che sto cercando, ma non sono sicuro se scipy ha già implementato questa funzionalità. Mi piacerebbe anche che funzioni per variabili non normali.

Sembra che lo Iman, Conover paper sia il metodo standard.

+0

È questo quello che stai cercando? http://stackoverflow.com/a/16025584/190597 – unutbu

+0

Funziona per variabili normali ... Ho altre distribuzioni. – MikeRand

+0

Sembra che il metodo raccomandato (Iman-Conover) usi una normale multi-variata per fare quello che sto cercando, quindi penso che il tuo commento sarà probabilmente un grande pezzo della soluzione finale (che è probabilmente qualcosa che farò costruire a mano). – MikeRand

risposta

8

Si desidera solo la correlazione attraverso una Copula gaussiana (*), quindi può essere calcolata in pochi passaggi con numpy e scipy.

  • creare variabili aleatorie multivariate con covarianza desiderato, numpy.random.multivariate_normal, e la creazione di un (nobs da k_variables) matrice

  • applica scipy.stats.norm.cdf per trasformare normale variabili casuali uniformi, per ogni colonna/variabile per ottenere uniforme marginale distribuzioni

  • applica dist.ppf trasformare margine uniforme la distribuzione desiderata, dove dist può essere una delle distribuzioni in scipy.stats

(*) copula gaussiana è solo una scelta e non è il migliore quando siamo interessati a un comportamento di coda, ma è il più facile da lavorare con ad esempio http://archive.wired.com/techbiz/it/magazine/17-03/wp_quant?currentPage=all

due riferimenti

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

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

(avrei potuto fare questo qualche tempo fa in python, ma non hanno script o funzioni in questo momento.)

+0

Conosci qualche soluzione simile, più efficiente in termini di memoria? Lo sto facendo con 'cov_matrix = toeplitz (rho ** arange (p))', ma mi imbatto in errori di memoria quando raggiungo dimensioni elevate. – MHankin

+0

Come posso ottenere distribuzioni marginali uniformi in python? – Ark

+0

@Ark Per ottenere distribuzioni marginali uniformi si salta l'ultimo passaggio. – user333700

3

Sembra che un metodo di campionamento basato sul rifiuto, come l'algoritmo di Metropolis-Hastings, sia ciò che si desidera. Scipy può implementare tali metodi con la sua funzione scipy.optimize.basinhopping.

I metodi di campionamento basati sul rifiuto consentono di prelevare campioni da una determinata distribuzione di probabilità. L'idea è di estrarre campioni casuali da un altro pdf "proposta" che è facile da campionare (come distribuzioni uniformi o gaussiane) e quindi utilizzare un test casuale per decidere se questo campione dalla distribuzione proposta deve essere "accettato" come rappresentante un campione della distribuzione desiderata.

I restanti trucchi saranno allora:

  1. figura il modulo della funzione articolare N-dimensionale densità di probabilità che ha marginali del modulo desiderato lungo ogni dimensione, ma con la matrice di correlazione che volere. Questo è facile da fare per la distribuzione gaussiana, dove la matrice di correlazione e il vettore medio desiderati sono tutto ciò che serve per definire la distribuzione. Se i tuoi marginali hanno un'espressione semplice, probabilmente puoi trovare questo pdf con un'algebra semplice ma noiosa. This cita molti altri che fanno quello di cui stai parlando, e sono certo che ce ne sono molti altri.

  2. Formulare una funzione per basinhopping in modo da ridurre al minimo tale che è accettata l'importo "minimo" per i campioni di questo pdf che è stato definito.

Dato il risultato di (1), (2) dovrebbe essere semplice.

Problemi correlati