Generating causally dependent random variables I'm trying to generate sets of causally connected random variables and started off doing this with a monte carlo approach.
The baseline is a 2-dimensional measured histogram from which I draw random values.
In my concrete examples these variables are acceleration $\bf{a}$ and velocity $\bf{v}$ - so obviously
$v_{i+1} = v_{i} + a_i * dt$
has to hold.
My current naive approach is:
I start with a some $v_0$.
Then I generate a random $a_0$ according to the measured probability of $\bf{a}$ for the value of $v_0$. Using this $a_0$ I can calculate $v_1$ and the whole procedure starts over again.
So when I check the generated accelerations $\bf{a}$ in bins of $\bf{v}$ everything's fine.
But I obviously this does not at all respect the marginal distribution of $\bf{v}$.
I'm kind of familiar with basic monte carlo methods, though lacking some theoretical background as you might guess.
I'd be fine if the two variables where just connected by some correlation matrix, but the causal connection between the two gives me headaches.
I didn't manage to find an example for this kind of problem somewhere - I might be googl'ing the wrong terms.
I'd be satisfied if somebody could point me to some literature/example or promising method to get a hold on this.
(Or tell me that's is not really possible given my inputs - that's what I'm guessing occasionally...)
EDIT:
The actual aim of this whole procedure:
I have a set of measurements $\bf{a}$ and $\bf{v}$, represented in a two-dimensional histogram $N(a,v)$. Given this input I'd like to generate sets of random $\bf{a_r}$ and $\bf{v_r}$ that reproduce the measured distribution.
 A: It seems that in order to reproduce the joint distribution $\rho(a,v)$, you should select new $a$ not only based on $v$, but based on the old $a$ also:
$a_{i+1} \sim \rho'(a_{i+1}|a_i, v_i)$
The question (to which I don't know the answer yet) is how to find $\rho'$ which produces $\rho$.
UPD:
You are to solve the following integral equation:
$$\rho(a, v) = \int da' \rho'\left(a|a', v-{a+a'\over 2}\Delta t\right) \rho(a', v-{a+a'\over 2}\Delta t)$$
Approximating the function $\rho$ with a histogram, you turn this to a system of linear equations:
$$\cases{
\rho(a, v) = \sum_{a'} \rho'\left(a|a', v-{a+a'\over 2}\Delta t\right) \rho(a', v-{a+a'\over 2}\Delta t) \\
\sum_a \rho'\left(a|a', v'\right) = 1}$$
This system is underdetermined. You may apply a smoothness penalty to obtain a solution.
A: Doesn't the gps data contain position $p$?  I would have thought that, not only is $v_{i+1}$ dependent upon $v_{i}$ and $a_{i}$ but $a_{i+1}$ would also be dependent upon $p_{i}$.  Consider:  in any road network there are bottlenecks, speed limits, signals, intersections, steep gradients, etc. that are geolocated.  So something like an ensemble (distribution) defined by:
$F_{a} = Pr ( A_{i+1} \le a_{i+1}\ |\ a_{i},v_{i},p_{i} )$
$v_{i+1} = v_{i} + a_{i}dt$  
For such an ensemble, the difficulty will lay in the nature of the data.  It is likely that the true population will be asymmetric, non-linear (piece-wise) and may not have defined moments.  These characteristics may not be evident within the sample you have at hand.
As @whuber has stated, the problem, ie exactly what you are seeking to produce, does not yet seem fully and clearly defined.  It is not clear as to whether you are interested in the ensemble or more so the individuals.
