Simulating the Poisson Point Process

The homogeneous Poisson point process is a counting process which satisfies {N(t),t>=0} where the probability of a random variable N being equal to n at time t is Poisson-distributed. The algorithm to simulate the Poisson point process can be simplified as follows (Chen et al. 2016).

The MATLAB script for simulating multivariate homogeneous Poisson point process is given below.

	function X = simpoisson(n, q, lambda, freq)
	% Simulation of poisson point processes
	%
	% Syntax:
	% X = simpoisson(n, q, lambda)
	%
	% Description:
	% It simulates a multivariate Poisson point process.
	%
	% Input Arguments:
	% n - the length of time series.
	% q - the dimension.
	% lambda - the parameter
	% freq - frequency
	%
	% Output Arguments:
	% X - a matrix of a multivariate Poisson point process.
	%
	% Examples:
	% X = simpoisson(2^15);
	% X = simpoisson(2^15, 2);
	% X = simpoisson(2^15, 10, 5);
	%
	% Copyright (c) 2017 Wonsang You(wsgyou@gmail.com)
	if nargin <2
	q = 1;
	lambda = 1;
	freq = 1;
	elseif nargin <3
	lambda = 1;
	freq = 1;
	elseif nargin <4
	freq = 1;
	end
	tmax = n*freq;
	tarr = zeros(1, q);
	% sums of exponentialy distributed interarrival times as matrix columns
	i = 1;
	while (min(tarr(i,:))<=tmax)
	tarr = [tarr; tarr(i, :)-log(rand(1, q))/lambda];
	i = i+1;
	end
	% cut off arrival times greater than tmax
	tarr(tarr>tmax) = tmax;
	t = 0:freq:tmax;
	X = zeros(n,q);
	for j = 1:q
	cumcnt = zeros(n,1);
	for i=1:n
	tarr_q = tarr(:,j);
	cumcnt(i)=length(find(tarr_q<t(i)));
	if i==1
	X(i,j) = cumcnt(i);
	else
	X(i,j) = cumcnt(i)-cumcnt(i-1);
	end
	end
	end

view raw simpoisson.m hosted with ❤ by GitHub

Run the following command on the MATLAB console to simulate a bivariate Poisson point process with λ=5 and length=128.

X = simpoisson(2^7,2,5); 
plot(X)

Then, you have two time series as displayed below.