# In a Poisson process measured with some efficiency, is the measured count still Poisson?

Situation:

Say I have a Poisson process, like radioactive decay, producing R particles per second. I measure with a detector. There is a probability P that a particle will be detected by the detector.

Things I think I know:

1. The inter-arrival time of the particle emission is exponentially distributed with parameters based on R.
2. The number of particles emitted before detection is given by a negative binomial based on P.
3. If a number N is sampled from (2), a single sample of inter-arrival time for detected particles can be given by the sum of N samples from (1). This sum can be obtained by sampling from a gamma distribution with parameters based on N and R.

My question:

If a single inter-arrival time can be calculated by sampling from a gamma based on N and R, how does the number of detector counts in an interval end up being Poisson again? (To be Poisson, the inter-arrival time for the detector must be exponential, not distributed according to some weird gamma thing.) Of course N is fluctuating, but I can't see how this works out.

However, I am almost completely sure detector counts are in fact Poisson distributed. Could somebody show me the math? Thanks for the help!

EDIT:

I found this paper: Fried, D. L. "Noise in photoemission current." Applied Optics 4.1 (1965): 79-80.

Which shows the result that a binomially selected poisson random variable is also Poisson with a rate given by PR. This confirms the comment by jbowman. Still, I would be interested in seeing the explanation of how my process of generating the inter-arrival time at the detector using the negative binomial and gamma distribution is incorrect. This is my major mental hiccup. Thank you.

EDIT 2:

I wrote this matlab script to test whether what I was trying with the gamma distribution worked. Turns out that somehow the gamma inter-arrival times generated with a geometrically distributed N are exponential and agree with the inter-arrival times suggested by Poisson(PR). (ia2 and ia3 are identically distributed). Any idea how this works out analytically? It was not intuitively obvious to me!

close all
n = 100000;
ia1 = exprnd(1,n,1); % create exponentially distributed inter-arrival times
t1 = cumsum(ia1); % running sum (the real experiment time)

mask = (rand(n,1) > 0.5); % flip a coin
t2 = t1(mask); % get only the events for which "the coin landed on heads"
ia2 = diff(t2); % calculate the inter-arrival times at the detector.

% plot the distributions
figure; hist(ia1,100); title('exponential inter-arrival times');
figure; hist(ia2,100); title('binomial sampled inter-arrival times');

%%
spacing = geornd(0.5,n,1) + 1; % how many events before we get heads
ia3 = gamrnd(spacing,ones(n,1)); % generate the interarrival times with gamma
figure; hist(ia3,100); title('geom/gamma inter-arrival times');

• #2 is actually not correct; if each particle has a probability $P$ of being detected, the distribution of detected particles per second is Poisson($RP$) (assuming detections are independent etc.) For #3, if $N$ is sampled from 2, then you don't have a single sample of inter-arrival times; you have a single observation of the sum of $N$ interarrival times, which is indeed distributed Gamma with shape parameter $N$. Consequently, the premise of your question ("If a single interarrival time...") isn't true. – jbowman May 29 '14 at 17:28
• I don't understand how you know the rate is Poisson(RP). Could you show me? That is the very heart of this question I think. In #2, I suppose that I have if I have a P chance to hit the detector, the number of particles emitted before hitting the detector is Geometrically distributed with a mean of 1/P. Thus, I can calculate can sample from this Geometric distribution to get N, then sum up N inter-arrival times to get a single inter-arrival time at the detector. Can you explain the flaw in this logic? I think your statement about the rate being Poisson(RP) is important. Thank you! – user487100 May 29 '14 at 18:14
• Are you somewhat familiar with moment generating / characteristic functions? I'd write it out using that approach, as it's simple, unless it's also unhelpful. – jbowman May 29 '14 at 18:21
• No I have not worked with moment generating functions. Do you have some idea of how to show that Poisson + some fixed probability of acceptance just scales the poisson rate? I am willing to learn the moment generating function approach based if you could show how this works out. – user487100 May 29 '14 at 18:43
• It'll be much later today (Pacific Std Time), I'm afraid; I can do it the straightforward way too, which will be less opaque. – jbowman May 29 '14 at 19:57

A quick non-technical argument might use Jackson networks. In your case total external arrivals is rate $R$, and there are no internal transitions (observed particles don't switch to the unobserved queue). The splitting proportion between the observed and unobserved nodes $p_{0i}$ is $P$, so the

$\lambda_{obs}=RP$

If you're looking for first principles, call $O(t)$ the observed counting process, and $N(t)\sim PP(r)$ the total counting process. Where each arrival in $N(t)$ gets logged in $O(t)$ with probability $p$. So that if for some $s$ we have $N(s)=n$ then $O(s)$ has a binomial($n,p$) distribution.

This approach uses probability generating functions:

$E[z^{O(t)}|N(t)=n]=\sum_{j=0}^{n}z^{j} {n \choose j}p^j(1-p)^{n-j}=(1-p+pz)^{n}$

Last equality by the binomial theorem. Then, unconditionally, since $N(t)\sim Poisson(rt)$:

$E[z^{O(t)}]=E[E[z^{O(t)}|N(t)=n]]=\sum_{n=0}^{\infty}(1-p+pz)^{n}\frac{rt^n}{n!}e^{-rt}=e^{-rt}e^{rt(1-p+pz)}=e^{rpt(z-1)}$

Which is the probability generating function of a Poisson($rpt$) random variable.