Serie storica simulata su R

Arado90
Qualcuno saprebbe linkarmi del materiale da cui "imparare" i comandi per simulare una serie storica con R?

Nel dettaglio, avrei bisogno di simulare una cosa del genere:

${(X_t\ |\ F_{t-1}:text{Poisson}(lambda_t); AAtinZZ) , (lambda_t=omega+alphaX_{t-1}):}$

(non mi viene la graffa grande su più righe :roll: )

In pratica la distribuzione dell'osservazione al tempo $t$ condizionata alle osservazioni passate è una Poisson con un parametro che dipende dall'osservazione precedente. Si tratta di un INARCH(1) che può essere riscritto come: $X_t=omega+alphaX_{t-1}+epsilon_t$ con $epsilon_t=X_t-lambda_t$

Io avrei bisogno di simulare un certo numero di osservazioni da un processo del genere utilizzando R.

Risposte
Arado90
Aggiornamento =D
Consultando alcuni testi/articoli ho trovato il seguente codice che permette di simulare un processo di Poisson non omogeneo:

t<-0
I<-0
T<-100
lambda<-1.5
r<-runif(1)
while(t<T){
    t=t-1/lambda*log(r)
    if(r<sin(t)/lambda)
    I<-c(I,t)
}


Il metodo consiste nello scegliere un $lambda>=lambda_t$, dove $lambda_t$ è la funzione d'intensità, che nell'esempio è la funzione seno mentre nel mio caso dovrebbe essere $omega+alpha*X_{t-1}$. E ogni $X_t$ dovrebbe essere un numero casuale estratto da una Poisson (esponenziale?) di parametro $lambda_t$ ($1/lambda_t$?).
Purtroppo ho fatto un po' di prove, ma non funziona :(
Tipo questa:

inarch<-function(omega,alfa,T,lambda){
  I<-0
  t<-0
  x<-NULL
  lambda<-NULL
  x[0]<-0
  for(t in 1:T){
    lambda[t]<-omega+alfa*x[t-1]
    x[t]<-rpois(1,lambda[t])
  }
  r<-runif(1)
  while(t<T){
    t=t-1/lambda*log(r)
    if(r<=lambda[t]/lambda)
    I<-c(I,t)
  }
}

hamming_burst
Ciao,
secondo me può essere utile che copi i risultati sbagliati oppure gli errori del debugger, che ottieni.
Provo a sposare la discussione in Informatica, vediamo se lì qualcuno che conosce R sa aiutarti.


lambda[t]/lambda

io cmq ti consiglierei di cambiare nome all'array lambda, non so come gestiche il typing R, di solito non è permesso l'overwriting nei linguaggi imperativi.

"Arado90":
E ogni $X_t$ dovrebbe essere un numero casuale estratto da una Poisson (esponenziale?) di parametro $lambda_t$ ($1/lambda_t$?).

che tipo di proprietà deve avere la distribuzione degli eventi?
EDIT: domanda idiota, lo hai scritto in altra maniera sopra :D

Arado90
io cmq ti consiglierei di cambiare nome all'array lambda, non so come gestiche il typing R, di solito non è permesso l'overwriting nei linguaggi imperativi.

C'hai ragione, correggiamo intanto questo!

inarch<-function(omega,alfa,T,lambda){
  I<-0
  t<-0
  x<-NULL
  lam<-NULL
  x[0]<-0
  for(t in 1:T){
    lam[t]<-omega+alfa*x[t-1]
    x[t]<-rpois(1,lam[t])
  }
  r<-runif(1)
  while(t<T){
    t=t-1/lambda*log(r)
    if(r<=lam[t]/lambda)
    I<-c(I,t)
  }
}

inarch(1,0.7,1000,1.5)
>Errore in rpois(1, lam[t]) : argomenti non validi



che tipo di proprietà deve avere la distribuzione degli eventi?

No, stavo pensando che $X_t$ è distribuita come una Poisson, però estraendo un solo numero casuale alla volta posso usare una esponenziale.

Il problema è che passando a rexp al posto di rpois il risultato finale è:
>There were 50 or more warnings (use warnings() to see the first 50)
xD

Tommyri
Guarda io di processi stocastici non ne so molto, tuttavia scrivere:

> x<-NULL
> x[0]<-0
> x
numeric(0)

genera un errore, mentre scrivere:
x<-NULL
> x[1]<-0
> x
[1] 0


credo faccia quanto voluto. Semplicemente i vettori in R cominciano dall'elemento 1 e non dall'elemento 0. Basta posticipare il ciclo di 1 e dovrebbe funzionare.

Arado90
Dannazione, grazie mille!! Mi era sfuggito questo particolare!
Adesso il programma quantomeno funziona =D

Resta da verificare la correttezza di quanto sto facendo, cioè se la serie di numeri che finalmente il programma restituisce è del tipo che serve a me.

Rispondi
Per rispondere a questa discussione devi prima effettuare il login.