Einstiegsschwierigkeiten, Daten aus Schleife nicht verlieren

[etwarda]

Hallo,

ich bin dabei gerade mein erstes kleines Projekt mit R zu starten doch stoße dabei immer wieder auf Schwierigkeiten.
Ich würde mich sehr freuen, wenn ihr mir ein paar Hinweise geben könntet.

Im Grunde möchte ich lediglich die folgende Gleichung (hier im Latex Code) in R implementieren:
c_{T}=-\frac{E((\sum_{s=0}^{T-1}a_{s})b_{T})}{E(b_{T})} mit T=1,...,50

Hierbei ist a ein zufällig generierter Vektor und b ergibt sich nach einer Bildungsvorschrift aus a.
In dem Beispiel ist c der 50-dimensionale Vektor, der sich aus den möglichen Werten von T ergibt. Der Erwartungswert E ergibt sich jeweils aus dem Mittelwert von 1000 Simulationen.

Ich bekomme es hin a und b aufzustellen. Allerdings bin ich unschlüssig, wie ich davon die 1000 Simulationen festhalten und in der Formel von c einbinden kann. Sollte ich dazu jeweils ein neuen Dataframe nutzen? Ich habe das Gefühl, dass ich mit den Tipps von viewtopic.php?f=20&t=30 auf dem richtigen Weg bin, doch bei jedem meiner Versuche werden mir nur Fehlermeldungen entgegen gebracht, weshalb ich gerade etwas ratlos bin.

Code: Alles auswählen

a<-numeric(50)
for (s in 1:50){
  a[s] <- rnorm(1, mean = 0, sd = 1)
}

b=numeric(50)
b[1]<-1-a[1]
for(s in 2:50){
  b[s]<-(1-a[s]*b[s-1])
}

c=numeric(50)
for(T in 1:50){
  c[T]<- ...
}

[jogo]

Hallo etwarda,

willkommen im Forum!

Code: Alles auswählen

a<-numeric(50)
for (s in 1:50){
  a[s] <- rnorm(1, mean = 0, sd = 1)
}

Das geht wesentlich kürzer mit

Code: Alles auswählen

T <- 50
a <- rnorm(T)

Bezüglich Simulation und Erwartungswert (bzw. arithmetischem Mittel) wird wohl die Funktion replicate() vorteilhaft eingesetzt werden können.

... ich bastel noch etwas weiter.

Nachtrag:
Die Formel ist doch äquivalent zu:
# c_{T} = -\frac{E( b_{T}(\sum_{s=0}^{T-1} a_{s}) )} {E(b_{T})}
# mit T=1,...,50
oder?

Dann könnte es etw so aussehen:

Code: Alles auswählen

sim <- function(T) {
  a <- rnorm(T)
  
  b <- numeric(T)
  b[1] <- 1-a[1]
  for(s in 2:T) b[s] <- 1-a[s]*b[s-1]
  
  
  c(Zaehler= b[T]*sum(a), bT=b[T])
}

T <- 50
f <- rowMeans(replicate(1000, sim(T)))
f

cT <- -f[1]/f[2]
cT

Gruß, Jörg