Proof for simulation of NHPP by thinning

Background:

I'm trying to show equivalency between the density function for a non-homogenous exponential process (NHEP?), (i.e. the arrival times of events generated by a non-homogenous Poisson process with time-varying rate parameter $$\lambda(t)$$) and the density function for arrival times generated by the "thinning" method for simulating a NHEP.

I know that the PDF of an NHEP is given by:

$$f_T(t) = \lambda(t) e^{-\int_0^t \lambda(t)dt}$$

I also know that an NHPP can be simulated by thinning, where we generate homogeneous exponential arrival times with rate parameter $$\lambda^* > \textrm{max}(\lambda(t))$$, and for each such arrival time, accept that arrival time with probability $$\frac{\lambda(t)}{\lambda^*}$$.

This process may be thought of as a combination of a standard exponential process and a Bernoulli random variable (call it $$p$$) that varies in time.

Question:

I'd like to somehow relate the thinning simulation method to the PDF of an NHPP; that is, show that the process produced by the thinning method is governed by the same underlying PDF.

Approach:

My idea is to show that the CDFs are the same. The CDF for the NHEP is given by:

$$\textrm{Pr}[T \leq t] = 1 - e^{-\int_0^t \lambda(t)dt}$$

So the question is, how to show that the CDF for the thinning simulation is equivalent.

To start, define:

• $$N^*_{T}$$ as the number of events generated by the "un-thinned," homogeneous exponential process over a given time-interval, $$T$$
• $$p_t$$ be the Bernoulli random variable at time $$t$$. If it $$p_t=0$$, the event at time $$t$$ (if present) is deleted. If $$p_t=1$$, the event at time $$t$$ (if present) is not deleted

I tried to write the CDF by considering the probability that no event occurred at each infinitesimally small time-step ($$\Delta t$$) between $$0$$ and $$T$$. At each such time-step, one of two conditions must be met:

1. No event is generated from the homogenous exponential process or
2. The Bernoulli random variable is zero

The probabilities of these two events are:

1. $$\textrm{Pr}[N^*_{\Delta t} = 0] = e^{-\lambda^* \Delta t}$$
2. $$\textrm{Pr}[p_t = 0] = 1 - \frac{\lambda(t)}{\lambda^*}$$

So, the probability that either of those events occur at a given infinitesimal time-step can be written as (skipping a few lines of algebra):

$$\textrm{Pr}[N^*_{\Delta t} = 0 \cup p_t=0] = 1 - \left(\frac{\lambda(t)}{\lambda^*}\left(1 - e^{-\lambda^* \Delta t}\right)\right)$$

From there, we can say that the probability of the first arrival time being greater than $$T$$ is the product of the probability of no events occurring during each of the infinitesimal time-steps between zero and $$T$$, i.e.

$$\textrm{Pr}[T > t] = \lim_{\Delta T \to 0} \prod_{i=0}^{\frac{t}{\Delta T}} \textrm{Pr}[N_{i} = 0 \cup p_{i\Delta T}=0]$$ $$\textrm{Pr}[T > t] = \lim_{\Delta T \to 0} \prod_{i=0}^{\frac{t}{\Delta T}} \left(1 - \left(\frac{\lambda(i\Delta T)}{\lambda^*}\left(1 - e^{-\lambda^* \Delta t}\right)\right)\right)$$

So then, the final CDF for the arrival times generated by the thinning process would be $$\textrm{Pr}[T \leq t] = 1-\textrm{Pr}[T > t] = 1-\lim_{\Delta T \to 0} \prod_{i=0}^{\frac{t}{\Delta T}} \left(1 - \left(\frac{\lambda(i\Delta T)}{\lambda^*}\left(1 - e^{-\lambda^* \Delta t}\right)\right)\right)$$

But, this is as far as I have gotten. I am at a loss as to how to show that this mess is equal to the NHEP CDF:

$$\textrm{Pr}[T \leq t] = 1 - e^{-\int_0^t \lambda(t)dt}$$

Or, perhaps this is the wrong approach altogether?

I think I've found an answer that at least satisfies me, so I thought I'd post it in case anyone else is interested in this question. It's a bit "hand-wavy" in parts, but I think it more or less makes sense. Picking up from where the question left off:

$$\textrm{Pr}[T \leq t] = 1-\lim_{\Delta T \to 0} \prod_{i=0}^{\frac{t}{\Delta T}} \left(1 - \left(\frac{\lambda(i\Delta T)}{\lambda^*}\left(1 - e^{-\lambda^* \Delta t}\right)\right)\right)$$

Now, this is one of the "hand-wavy" bits. Based on the Taylor series for $$e^x$$, we can say:

$$e^{-\lambda^* \Delta t} = 1 + (-\lambda^* \Delta t) + \frac{(-\lambda^* \Delta t)^2}{2!} + \frac{(-\lambda^* \Delta t)^3}{3!} + ...$$

But, we also know that we're dealing with the limit as $$\Delta t \to 0$$. So, I'm thinking we can neglect all but the first order terms.

$$\lim_{\Delta t \to 0} e^{-\lambda^* \Delta t} = 1 -\lambda^* \Delta t$$

Using this expression, we can rewrite the probability expression as:

$$\textrm{Pr}[T \leq t] = 1-\lim_{\Delta T \to 0} \prod_{i=0}^{\frac{t}{\Delta T}} \left(1 - \left(\frac{\lambda(i\Delta T)}{\lambda^*}\left(1 - 1 + \lambda^* \Delta t\right)\right)\right)\\ = 1-\lim_{\Delta T \to 0} \prod_{i=0}^{\frac{t}{\Delta T}} \left(1 - \lambda(i\Delta T)\Delta t\right)\\ = 1-\exp\left(\lim_{\Delta T \to 0} \ln \prod_{i=0}^{\frac{t}{\Delta T}} \left(1 - \lambda(i\Delta T)\Delta t\right)\right)\\ = 1-\exp\left(\lim_{\Delta T \to 0} \sum_{i=0}^{\frac{t}{\Delta T}} \ln \left(1 - \lambda(i\Delta T)\Delta t\right)\right)$$

Now, we can do the reverse approximation. Specifically, $$\ln(1+x) = x + O(x^2)$$ as $$x\to 0$$, so we can say:

$$\textrm{Pr}[T \leq t] = 1-\exp\left(\lim_{\Delta T \to 0} \sum_{i=0}^{\frac{t}{\Delta T}} \lambda(i\Delta T)\Delta t\right)$$

But now the term in the integral is just the Riemann sum for the integral, i.e.

$$\lim_{\Delta T \to 0} \sum_{i=0}^{\frac{t}{\Delta T}} \lambda(i\Delta T)\Delta t = \int_0^T \lambda(t) dt$$

So: $$\textrm{Pr}[T \leq t] = 1-e^{\int_0^T \lambda(t) dt}$$