I have a set $X$ of $10^6$ elements and a time series of probability distributions $\mu_1,\mu_2,\ldots$ on $X$. I want to estimate the expected value of a function $f$ over each $\mu_t$. It is easy to evaluate $\mu_t$, but it is expensive to evaluate $f$, so I want to use a set of $10^3$ (possibly weighted) samples from $\mu_t$.
So far, so vanilla. Here's the twist: Each $\mu_{t+1}$ is highly correlated with the preceding $\mu_t$, so I think there is an opportunity to reuse samples. The question is: how?
Specifically, I want to generate a sequence of sample sets $S_1\subseteq S_2\subseteq\cdots \subsetneq X$ and a sequence of weight functions $w_t:S_t\to\mathbb R$ such that the estimate of $\sum_{x\in X}f(x)\mu_t(x)$ is given by $\sum_{x\in S_t}f(x)w_t(x)$.
Requirements:
- The algorithm should call for $<<10^3$ new samples on every time step. Ideally it should dynamically adapt to the degree of correlation between distributions.
- The samples and weights at time $t$ must not depend on future times.
- The algorithm must still "work" in the case that an element $x$ has $\mu_1(x)=0$ but $\mu_t(x)$ eventually becomes large.
- I don't mind if the estimates are correlated with each other to a greater degree than estimates from independent samples would be correlated. What I want to ensure is that the estimate for a given time is not biased by the distributions for preceding times.
For example, I could let $S_1$ be a sample of $N$ elements chosen from $\mu_1$ and initialize $w_1(x)=1/N$. Then, for each subsequent time $t$, reuse $S_t=S_1$ and just reweight $w_t(x)=\mu_t(x)/\mu_1(x)$ (and then normalize $w_t$). This algorithm is perfectly efficient, in that the sample set never grows and $f$ never has to be evaluated after the first estimate. It passes requirements (1) and (2) with flying colors. However, I believe it fails (3) and (4). I'm not sure how to precisely express the way in which the algorithm fails those requirements, so feel free to enlighten me on that as well!
My intuition is that a solution will look something like this: Let $S_{t+1}\setminus S_t$ be a set of samples from $\max(0,\mu_{t+1}-\mu_t)$ with constant weights, and also import $S_t$ into $S_{t+1}$ with weights adjusted as in the previous paragraph. But that seems too ad-hoc to be correct, and it doesn't say how to choose the number of new samples that are needed to maintain the quality (variance?) of the estimate.
Perhaps if I decrease the old weights by a global factor of $\min_{x\in X}(\mu_{t+1}/\mu_t)$? But in that case, I'm concerned that the minimum will tend to be too small because some $x$ has an unlucky decrease in probability, so each update will almost have to start from scratch. Would it be good enough to decrease by $\min_{x\in S_t}(\mu_{t+1}/\mu_t)$?
(I added the resampling tag because my idea is reminiscent of the weight-adjusting step in a particle filter, which is followed by resampling and then, periodically, refreshing the sample set with new samples. I am not sure if there is a precise mapping between my problem and particle filters, but feel free to explore that connection in an answer.)
