Question About Sampling from Conditional Distribution Created from Approximated Distribution

Question

First, set the data and prepare the conditional distribution.


pxy is joint distribution of Random Variable x andy.
py is marginal distribution for variable y of pxy
pcxy is conditional distribution(distribution when x==0)

data = RandomFunction[WienerProcess[0, 1], {0, 1, 0.01}, 1000];
data // Normal // #[[All, 2]] & /@ # & // Set[sample, #] &;

Table[Transpose@{sample[[All, i]], sample[[All, i + 1]]}, {i, 1, 
      99}] // Flatten[#, 1] & // KernelMixtureDistribution // 
  Set[pxy, #] &;
px = MarginalDistribution[pxy, 1];
pcxy = PDF[pxy, {0, y}]/PDF[px, 0];

What I want to do is sampling from pcxy quickly,or fast as long as possible.

Plot[pcxy, {y, -1, 1}]

What I've tried is the following(this didn't end at all):

ProbabilityDistribution[pcxy, {y, -Infinity, Infinity}]
RandomVariate[%, 3] // AbsoluteTiming

Answer 1

The reason it's slow is that the PDF of KernelMixture distributions is given symbolically. This can lead to very large expressions if you're using a lot of points, which would then need to be integrated for RandomVariate to work.

Using an interpolated PDF is a decent alternative

data = RandomFunction[WienerProcess[0, 1], {0, 1, 0.01}, 1000];

sample = data["ValueList"];

pxy = KernelMixtureDistribution[
    Flatten[Transpose /@ Partition[Transpose[sample[[;;, ;; 100]]], 2, 1], 1]
];

px = MarginalDistribution[pxy, 1];

pcxypdf = Interpolation[
    Table[{y, PDF[pxy, {0.0, y}] / PDF[px, 0.0]}, {y, -1, 1, 0.01}]
];

pcxy = ProbabilityDistribution[pcxypdf[y], {y, -1, 1}];

RandomVariate[pcxy, 1000] // AbsoluteTiming // First
(* 0.062171 *)

Show[
    Plot[PDF[pxy, {0, y}] / PDF[px, 0], {y, -1, 1}, PlotStyle -> Red],
    SmoothHistogram[RandomVariate[pcxy, 1000], Automatic, "PDF"],
    ImageSize -> Large
]

Answer 2

One expects the conditional distribution to look very much like a normal distribution and indeed your plot of the estimated density shows that. Therefore, you could sample from a normal distribution with the same mean and standard deviation as the estimated conditional distribution.

data = RandomFunction[WienerProcess[0, 1], {0, 1, 0.01}, 1000];
data // Normal // #[[All, 2]] & /@ # & // Set[sample, #] &;

Table[Transpose@{sample[[All, i]], sample[[All, i + 1]]}, {i, 1, 
      99}] //
    Flatten[#, 1] & // KernelMixtureDistribution // Set[pxy, #] &;
px = MarginalDistribution[pxy, 1];
pcxy = PDF[pxy, {0, y}]/PDF[px, 0];

ymean = NIntegrate[y pcxy, {y, -1, 1}]
ysd = NIntegrate[(y - ymean)^2 pcxy, {y, -1, 1}]^0.5
Plot[{pcxy, PDF[NormalDistribution[ymean, ysd], y]}, {y, -1, 1},
 PlotStyle -> {{Thickness[0.02], LightGray}, {Thin, Red}}]

Future samples could then be found with

yFuture = RandomVariate[NormalDistribution[ymean, ysd], 100];