I have the following probability model:
$(X_k|\text{ PastHistory}_{k-1}, \theta_0, \theta_1, \theta_2) \sim (\pi\cdot N(\theta_1+\theta_0\cdot X_{k-1},1)+(1-\pi)\cdot N(\theta_2+\theta_0\cdot X_{k-1},1))$
and $X_k=0.5X_{k-1}+\epsilon_k$ with $\epsilon_k \sim N(0,1)$
I've tried something like:
Sample[n_] :=
Module[{lst},
lst = {RandomVariate[NormalDistribution[0, 1/(1 - 0.5^2)]]};
element = 1;
While[element <= n,
AppendTo[lst,
RandomVariate[
TransformedDistribution[
p*X1 + (1 - p) X2, {
X1 \[Distributed]
NormalDistribution[θ1 + θ0*lst[[element]], 1],
X2 \[Distributed]
NormalDistribution[θ2 + θ0*lst[[element]], 1]
}]]];
element = element + 1;
];
lst
];
Is this correct?
AppendTois very slow. See also the Q&A something faster than AppendTo in a loop – Jacob Akkerboom Mar 03 '14 at 16:46TransformedDistribution[p*X1 + (1 - p) X2, {X1 \[Distributed] NormalDistribution[θ1 + θ0*el, 1], X2 \[Distributed] NormalDistribution[θ2 + θ0*el, 1]}]simplifies toNormalDistribution[el θ0+p (θ1-θ2)+θ2,Sqrt[1+2 (-1+p) p]]. It is better to put the simplified version in the code, otherwise this simplification will be done at each iteration. – Jacob Akkerboom Mar 03 '14 at 17:01θ0 = 1/2for instance. You seem to only use the first equation for $X_k$ in your code. – Jacob Akkerboom Mar 03 '14 at 17:13