FAQFAQ   SearchSearch   MemberlistMemberlist   UsergroupsUsergroups 
 ProfileProfile   PreferencesPreferences   Log in to check your private messagesLog in to check your private messages   Log inLog in 
Forum index » Science and Technology » Math » Probability
Poisson process: analysis and simulation disagree
Post new topic   Reply to topic Page 1 of 1 [15 Posts] View previous topic :: View next topic
Author Message
ijs@iitis.gliwice.pl
science forum beginner


Joined: 29 Apr 2006
Posts: 10

PostPosted: Sat Apr 29, 2006 1:26 am    Post subject: Poisson process: analysis and simulation disagree Reply with quote

Hi,

I have a task which involves a Poisson process, and I want to
calculate certain probabilities. However, my simulation and analysis
return different results, and so I must be making some mistake in this
rather simple task. I would appreciate it if someone could write
where I'm making my mistake. First I will describe my analytical
solution, then I will write about my simulation.


PROBLEM STATEMENT
-----------------

A doctor can accept one patient per hour. Patients arrive at the
doctor's office according to the Poisson distribution with some given
intensity u. So the probability that there are k patients is:

f(k, u) = e^(-u) * u^k / k!

But the doctor can accept only one patient chosen at random, and the
rest has to go home. The rejected patients don't wait for the next
hour. Therefore the probability of being accepted when there are k
patients equals 1/k.

My question is this: what's the probability that I will be accepted
when I go to see this doctor?


ANALYSIS
--------

Event A: a patient is accepted
Event B: at least one patient arrives
Event B_k: k patients arrive

B = B_1 or B_2 or ... or B_n

where "or" is the set sum operator, and n is so big, that f(n, u) is
small enough (for u = 1.0, n = 10 is enough).

We want to get the probability P(A|B), that is the probability that a
patient is accepted provided that at least one patient arrived.

P(A|B) = P(A and B) / P(B)

where "and" is the set product operator.

P(B) = sum of f(k, u) for k = 1 to n

P(A and B) = P(A and (B_1 or B_2 or ... or B_n))

Because B_k and pair-wise disjoint, i.e. (B_i and B_j = 0, for i !=
j), then:

P(A and (B_1 or B_2 or ... or B_n)) =
P(A and B_1) + ... + P(A and B_n)

Moreover:

P(A and B_k) = P(A|B_k) * P(B_k)

The above we can compute having:

P(A|B_k) = 1 / k
P(B_k) = f(k, u)

Therefore the complete expression for the probability P(A|B) that we
look for is this:

sum of ((1 / k) * f(k, u)) for k = 1 to n
P(A|B) = -----------------------------------------
sum of f(k, u) for k = 1 to n

That's it for analysis.


PROGRAM FOR ANALYSIS
--------------------

In this section and the next one I will give two programs. They both
require the GNU Scientific Library (GSL). I compiled and ran these
programs on Linux 2.6.10, GCC compiler v. 3.4.4, and GSL 1.5. My
Makefile (works with GNU Make) for building both programs is this:

all: ana sim

LDFLAGS = -l gsl -l gslcblas


My program for analysis is this:

#include <stdio.h>
#include <gsl/gsl_rng.h>
#include <gsl/gsl_randist.h>

int
main (void)
{
const gsl_rng_type * T;
gsl_rng * r;
int k, n = 10;
double u = 1.0;
double nominator = 0;
double denominator = 0;

gsl_rng_env_setup();
T = gsl_rng_default;
r = gsl_rng_alloc (T);

for (k = 1; k < n; k++)
{
nominator += (1.0 / (double)k) * gsl_ran_poisson_pdf(k, u);
denominator += gsl_ran_poisson_pdf(k, u);
}

printf("u = %f\n", u);
printf("P(A|B) = %f\n", nominator / denominator);
return 0;
}

In this program the function f(k, u) is gsl_ran_poisson_pdf(k, u).


SIMULATION
----------

The simulation program is this:

#include <stdio.h>
#include <gsl/gsl_rng.h>
#include <gsl/gsl_randist.h>

int
main (void)
{
const gsl_rng_type * T;
gsl_rng * r;
int i, k;
int arrived = 0;
int accepted = 0;
int max_i = 10000000;
double u = 1.0;

gsl_rng_env_setup();
T = gsl_rng_default;
r = gsl_rng_alloc (T);

// In one iteration of this loop we simulate one hour.
for (i = 0; i < max_i; i++)
{
k = gsl_ran_poisson(r, u);

if (k)
{
arrived += k;
accepted += 1;
}
}

printf("u = %f\n", u);
printf("P(A|B) = %f\n", (double) accepted / (double) arrived);

return 0;
}

The function gsl_ran_poisson(r, u) returns a random integer according
to the Poisson distribution with the intensity u.


COMPARISON OF SIMULATION AND ANALYSIS
-------------------------------------

When you run the analysis, you'll get:

u = 1.000000
P(A|B) = 0.766988

while with the simulation:

u = 1.000000
P(A|B) = 0.632072


The two reported values of P(A|B) should be the same, and I'm unable
to explain why they are not. I would appreciate any advice.

Thanks for reading.


Best,
Irek
Back to top
Dan Akers
science forum addict


Joined: 19 Jul 2005
Posts: 56

PostPosted: Sat Apr 29, 2006 2:33 am    Post subject: Re: Poisson process: analysis and simulation disagree Reply with quote

You wrote:
"I have a task which involves a Poisson process, and I want to calculate
certain probabilities. However, my simulation and analysis return
different results, and so I must be making some mistake in this rather
simple task. I would appreciate it if someone could write where I'm
making my mistake. First I will describe my analytical solution, then I
will write about my simulation.
PROBLEM STATEMENT
-----------------
A doctor can accept one patient per hour. Patients arrive at the
doctor's office according to the Poisson distribution with some given
intensity u. So the probability that there are k patients is:
f(k, u) = e^(-u) * u^k / k!
But the doctor can accept only one patient chosen at random, and the
rest has to go home. The rejected patients don't wait for the next hour.
Therefore the probability of being accepted when there are k patients
equals 1/k.
My question is this: what's the probability that I will be accepted when
I go to see this doctor?<snip>
______________________________________
Re;
Here's my stab at it: I tend to agree with your simulation and not your
analysis. Given that u (mean patient arrival rate) is 1/hr and
presumably, it takes one hour for the doctor to see a patient. By
definition then, over a long period of time (t), there will be (t)(u)
patients who attempt to see the doctor. However, Poisson tells us that
the probability of 0 patients showing up for any given t-unit (hourly in
this case) interval is 0.3678. It follows then, that 1-0.3678 or 63.22%
of the patients who go must get to see the doctor. This also implies
that the doctor is idle for 36.78% of the time because of such a
schedule regimen.

Dan Akers
Back to top
ijs@iitis.gliwice.pl
science forum beginner


Joined: 29 Apr 2006
Posts: 10

PostPosted: Sat Apr 29, 2006 3:54 am    Post subject: Re: Poisson process: analysis and simulation disagree Reply with quote

Dan,

Thanks for your post.

Quote:
Here's my stab at it: I tend to agree with your simulation and not
your analysis.

I agree with you: the simulation is less complicated, and I also trust
it more.

Quote:
Given that u (mean patient arrival rate) is 1/hr and presumably,
it takes one hour for the doctor to see a patient. By definition
then, over a long period of time (t), there will be (t)(u) patients
who attempt to see the doctor. However, Poisson tells us that
the probability of 0 patients showing up for any given t-unit
(hourly in this case) interval is 0.3678. It follows then, that
1-0.3678 or 63.22% of the patients who go must get to see
the doctor. This also implies that the doctor is idle for 36.78%
of the time because of such a schedule regimen.

Agreed. I wonder where in my analysis I'm making an error.


Best,
Irek
Back to top
Dan Akers
science forum addict


Joined: 19 Jul 2005
Posts: 56

PostPosted: Sat Apr 29, 2006 2:35 pm    Post subject: Re: Poisson process: analysis and simulation disagree Reply with quote

Irek wrote:
"Agreed. I wonder where in my analysis I'm making an error."
__________________________________
Re:
I think it may be the conditional statement. I don't think you need it.
In fact, I was considering that my answer of 63.22% was too high. I say
this, because I know that over a great number of time periods, say 1000
hrs, that Poisson estimates that 36.68% of them will have no patients
show up, given the mean arrival rate of 1/hr. That means that out of
1000 hrs, 369 or so, have no patient arrivals. Since the doctor also
can only see patients at the rate of 1/hr as well, it follows then that
for about 369 out of every 1000 hrs, he is idle. So, as an estimate,
the doctor can only see 631 or so patients per 1000 hrs. And if you are
one of the 1000 patients that arrives during that 1000 hr time period,
you have about a 63.2% chance of seeing the doctor and not being turned
away.

Dan Akers
Back to top
Dan Akers
science forum addict


Joined: 19 Jul 2005
Posts: 56

PostPosted: Sat Apr 29, 2006 2:58 pm    Post subject: Re: Poisson process: analysis and simulation disagree Reply with quote

Irek wrote:
"But the doctor can accept only one patient chosen at random, and the
rest has to go home. The rejected patients don't wait for the next hour.
Therefore the probability of being accepted when there are k patients
equals 1/k."
______________________________________
Re;
This is the key I think. My "analysis" assumes that there is no queuing
of patients; no patients waiting; therefor the doctor does not "choose"
the next patient. Hence all that idle time for the doctor. As you said
in the problem statement, patients who arrive are turned away
immediately if the doctor already has a patient. Thus, there is only
one random variable in the process and that is the arrival of patients,
no doctor choice. However, your answer of 0.766 from you analysis may
indeed be correct for queued patients, since allowing the patients to
queue would certainly cut the doctor's idle time. I'm still thinking
about that...

Dan Akers
Back to top
ijs@iitis.gliwice.pl
science forum beginner


Joined: 29 Apr 2006
Posts: 10

PostPosted: Sat Apr 29, 2006 5:06 pm    Post subject: Re: Poisson process: analysis and simulation disagree Reply with quote

This version of the simulation returns both results: one that agrees
with the analysis, the other that doesn't. The latter result seems to
be right. I expected both results to be equal. I still don't know
where in the analysis I'm making a mistake.

Best,
Irek

*************************************************************

#include <stdio.h>
#include <gsl/gsl_rng.h>
#include <gsl/gsl_randist.h>

#define N 100

int arrived[N];

int
main (void)
{
const gsl_rng_type * T;
gsl_rng * r;
int i, k;
int max_i = 10000000;
double u = 1.0;
double result = 0;
double arr_sum;

int num_arrived = 0;
int num_routed = 0;

gsl_rng_env_setup();
T = gsl_rng_default;
r = gsl_rng_alloc (T);

for (i = 0; i < max_i; i++)
{
k = gsl_ran_poisson(r, u);

if (k && k < N)
{
arrived[k]++;
num_arrived += k;
num_routed++;
}
}

for (i = 1; i < N; i++)
{
arr_sum += arrived[i];
result += 1.0 / i * (double) arrived[i] / max_i;
}

arr_sum /= max_i;

printf("u = %f\n", u);
printf("P(A|B) = %f\n", result / arr_sum);
printf("P(A|B) = %f\n", (double) num_routed / num_arrived);

return 0;
}
Back to top
Mike1170
science forum addict


Joined: 17 Sep 2005
Posts: 74

PostPosted: Sun Apr 30, 2006 12:45 am    Post subject: Re: Poisson process: analysis and simulation disagree Reply with quote

Can't follow your analysis and don't see why you define events A and B as
you do. You only need to focus on the hour in which YOU arrive in which
case the doctor will ALWAYS find at least one patient to choose from--YOU.
What is the point of defining Event A as "a patient is selected" and B as
"at least one patient arrives"?
P{A} = 1 and P{B} = 1.

Like Dan A, I agree with your simulation (not your analysis) and think of
the situation as follows:

Doctor opens his door once per hour and randomly selects one of the patients
waiting.
Patients not selected disappear instantly and new patients begin to arrive
and accumulate over the next hour until the doctor opens his door again to
make another selection. In this context, the Poisson distribution is
sensible for the number of patients (besides yourself) who have accumulated
in the waiting room over the hour, in that it allows the inter-arrival times
of patients to the waiting room to be exponentially distributed (average
rate of arrival = 1/u). (The exponential distribution is commonly used to
describe the time between successive arrivals of whatever and the Poisson
the number of such arrivals in a fixed period of time).

You, in effect, join a group whose size is (already) a Poisson variable.
Thus, if the "group" you join happens to be of size 0, your probability of
being selected is 100%. If the group you join is of size 1, your probability
of being selected is 50% and so on.

Your probability of being selected is exp(-u) times the summation of
(u^k)/(k+1)! from k=0 to infinity.
I didn't try to simply the series but I think it can be.
Anyway, a quick simulation produced the following (value of p agrees with
your simulated value when u = 1)
u= 0.01, p=0.9950;
u= 0.1, p=0.9516;
u= 0.5 p=0.7869;
u=1.0, p=0.6321;
u=2.0, p=0.4323;
u=3.0, p=0.3167;
u=10, p=0.1.
p seems to approach 1/u as u gets larger. Not sure why it doesn't approach
1/(u+1).



<ijs@iitis.gliwice.pl> wrote
Quote:
Hi,

I have a task which involves a Poisson process, and I want to
calculate certain probabilities. However, my simulation and analysis
return different results, and so I must be making some mistake in this
rather simple task. I would appreciate it if someone could write
where I'm making my mistake. First I will describe my analytical
solution, then I will write about my simulation.


PROBLEM STATEMENT
-----------------

A doctor can accept one patient per hour. Patients arrive at the
doctor's office according to the Poisson distribution with some given
intensity u. So the probability that there are k patients is:

f(k, u) = e^(-u) * u^k / k!

But the doctor can accept only one patient chosen at random, and the
rest has to go home. The rejected patients don't wait for the next
hour. Therefore the probability of being accepted when there are k
patients equals 1/k.

My question is this: what's the probability that I will be accepted
when I go to see this doctor?


ANALYSIS
--------

Event A: a patient is accepted
Event B: at least one patient arrives
Event B_k: k patients arrive

B = B_1 or B_2 or ... or B_n

where "or" is the set sum operator, and n is so big, that f(n, u) is
small enough (for u = 1.0, n = 10 is enough).

We want to get the probability P(A|B), that is the probability that a
patient is accepted provided that at least one patient arrived.

P(A|B) = P(A and B) / P(B)

where "and" is the set product operator.

P(B) = sum of f(k, u) for k = 1 to n

P(A and B) = P(A and (B_1 or B_2 or ... or B_n))

Because B_k and pair-wise disjoint, i.e. (B_i and B_j = 0, for i !=
j), then:

P(A and (B_1 or B_2 or ... or B_n)) =
P(A and B_1) + ... + P(A and B_n)

Moreover:

P(A and B_k) = P(A|B_k) * P(B_k)

The above we can compute having:

P(A|B_k) = 1 / k
P(B_k) = f(k, u)

Therefore the complete expression for the probability P(A|B) that we
look for is this:

sum of ((1 / k) * f(k, u)) for k = 1 to n
P(A|B) = -----------------------------------------
sum of f(k, u) for k = 1 to n

That's it for analysis.


PROGRAM FOR ANALYSIS
--------------------

In this section and the next one I will give two programs. They both
require the GNU Scientific Library (GSL). I compiled and ran these
programs on Linux 2.6.10, GCC compiler v. 3.4.4, and GSL 1.5. My
Makefile (works with GNU Make) for building both programs is this:

all: ana sim

LDFLAGS = -l gsl -l gslcblas


My program for analysis is this:

#include <stdio.h
#include <gsl/gsl_rng.h
#include <gsl/gsl_randist.h

int
main (void)
{
const gsl_rng_type * T;
gsl_rng * r;
int k, n = 10;
double u = 1.0;
double nominator = 0;
double denominator = 0;

gsl_rng_env_setup();
T = gsl_rng_default;
r = gsl_rng_alloc (T);

for (k = 1; k < n; k++)
{
nominator += (1.0 / (double)k) * gsl_ran_poisson_pdf(k, u);
denominator += gsl_ran_poisson_pdf(k, u);
}

printf("u = %f\n", u);
printf("P(A|B) = %f\n", nominator / denominator);
return 0;
}

In this program the function f(k, u) is gsl_ran_poisson_pdf(k, u).


SIMULATION
----------

The simulation program is this:

#include <stdio.h
#include <gsl/gsl_rng.h
#include <gsl/gsl_randist.h

int
main (void)
{
const gsl_rng_type * T;
gsl_rng * r;
int i, k;
int arrived = 0;
int accepted = 0;
int max_i = 10000000;
double u = 1.0;

gsl_rng_env_setup();
T = gsl_rng_default;
r = gsl_rng_alloc (T);

// In one iteration of this loop we simulate one hour.
for (i = 0; i < max_i; i++)
{
k = gsl_ran_poisson(r, u);

if (k)
{
arrived += k;
accepted += 1;
}
}

printf("u = %f\n", u);
printf("P(A|B) = %f\n", (double) accepted / (double) arrived);

return 0;
}

The function gsl_ran_poisson(r, u) returns a random integer according
to the Poisson distribution with the intensity u.


COMPARISON OF SIMULATION AND ANALYSIS
-------------------------------------

When you run the analysis, you'll get:

u = 1.000000
P(A|B) = 0.766988

while with the simulation:

u = 1.000000
P(A|B) = 0.632072


The two reported values of P(A|B) should be the same, and I'm unable
to explain why they are not. I would appreciate any advice.

Thanks for reading.


Best,
Irek
Back to top
Mike1170
science forum addict


Joined: 17 Sep 2005
Posts: 74

PostPosted: Sun Apr 30, 2006 12:59 am    Post subject: Re: Poisson process: analysis and simulation disagree Reply with quote

Further to comments below, series simplifies to [1-exp(-u)] / u (=
probability of being selected by doctor).

"Mike" <MDC@nomail.com> wrote in message news:6RT4g.40$pP7.32@fe11.lga...
Quote:
Can't follow your analysis and don't see why you define events A and B as
you do. You only need to focus on the hour in which YOU arrive in which
case the doctor will ALWAYS find at least one patient to choose from--YOU.
What is the point of defining Event A as "a patient is selected" and B as
"at least one patient arrives"?
P{A} = 1 and P{B} = 1.

Like Dan A, I agree with your simulation (not your analysis) and think of
the situation as follows:

Doctor opens his door once per hour and randomly selects one of the
patients waiting.
Patients not selected disappear instantly and new patients begin to arrive
and accumulate over the next hour until the doctor opens his door again to
make another selection. In this context, the Poisson distribution is
sensible for the number of patients (besides yourself) who have
accumulated in the waiting room over the hour, in that it allows the
inter-arrival times of patients to the waiting room to be exponentially
distributed (average rate of arrival = 1/u). (The exponential
distribution is commonly used to describe the time between successive
arrivals of whatever and the Poisson the number of such arrivals in a
fixed period of time).

You, in effect, join a group whose size is (already) a Poisson variable.
Thus, if the "group" you join happens to be of size 0, your probability of
being selected is 100%. If the group you join is of size 1, your
probability of being selected is 50% and so on.

Your probability of being selected is exp(-u) times the summation of
(u^k)/(k+1)! from k=0 to infinity.
I didn't try to simply the series but I think it can be.
Anyway, a quick simulation produced the following (value of p agrees with
your simulated value when u = 1)
u= 0.01, p=0.9950;
u= 0.1, p=0.9516;
u= 0.5 p=0.7869;
u=1.0, p=0.6321;
u=2.0, p=0.4323;
u=3.0, p=0.3167;
u=10, p=0.1.
p seems to approach 1/u as u gets larger. Not sure why it doesn't
approach 1/(u+1).



ijs@iitis.gliwice.pl> wrote
Hi,

I have a task which involves a Poisson process, and I want to
calculate certain probabilities. However, my simulation and analysis
return different results, and so I must be making some mistake in this
rather simple task. I would appreciate it if someone could write
where I'm making my mistake. First I will describe my analytical
solution, then I will write about my simulation.


PROBLEM STATEMENT
-----------------

A doctor can accept one patient per hour. Patients arrive at the
doctor's office according to the Poisson distribution with some given
intensity u. So the probability that there are k patients is:

f(k, u) = e^(-u) * u^k / k!

But the doctor can accept only one patient chosen at random, and the
rest has to go home. The rejected patients don't wait for the next
hour. Therefore the probability of being accepted when there are k
patients equals 1/k.

My question is this: what's the probability that I will be accepted
when I go to see this doctor?


ANALYSIS
--------

Event A: a patient is accepted
Event B: at least one patient arrives
Event B_k: k patients arrive

B = B_1 or B_2 or ... or B_n

where "or" is the set sum operator, and n is so big, that f(n, u) is
small enough (for u = 1.0, n = 10 is enough).

We want to get the probability P(A|B), that is the probability that a
patient is accepted provided that at least one patient arrived.

P(A|B) = P(A and B) / P(B)

where "and" is the set product operator.

P(B) = sum of f(k, u) for k = 1 to n

P(A and B) = P(A and (B_1 or B_2 or ... or B_n))

Because B_k and pair-wise disjoint, i.e. (B_i and B_j = 0, for i !=
j), then:

P(A and (B_1 or B_2 or ... or B_n)) =
P(A and B_1) + ... + P(A and B_n)

Moreover:

P(A and B_k) = P(A|B_k) * P(B_k)

The above we can compute having:

P(A|B_k) = 1 / k
P(B_k) = f(k, u)

Therefore the complete expression for the probability P(A|B) that we
look for is this:

sum of ((1 / k) * f(k, u)) for k = 1 to n
P(A|B) = -----------------------------------------
sum of f(k, u) for k = 1 to n

That's it for analysis.


PROGRAM FOR ANALYSIS
--------------------

In this section and the next one I will give two programs. They both
require the GNU Scientific Library (GSL). I compiled and ran these
programs on Linux 2.6.10, GCC compiler v. 3.4.4, and GSL 1.5. My
Makefile (works with GNU Make) for building both programs is this:

all: ana sim

LDFLAGS = -l gsl -l gslcblas


My program for analysis is this:

#include <stdio.h
#include <gsl/gsl_rng.h
#include <gsl/gsl_randist.h

int
main (void)
{
const gsl_rng_type * T;
gsl_rng * r;
int k, n = 10;
double u = 1.0;
double nominator = 0;
double denominator = 0;

gsl_rng_env_setup();
T = gsl_rng_default;
r = gsl_rng_alloc (T);

for (k = 1; k < n; k++)
{
nominator += (1.0 / (double)k) * gsl_ran_poisson_pdf(k, u);
denominator += gsl_ran_poisson_pdf(k, u);
}

printf("u = %f\n", u);
printf("P(A|B) = %f\n", nominator / denominator);
return 0;
}

In this program the function f(k, u) is gsl_ran_poisson_pdf(k, u).


SIMULATION
----------

The simulation program is this:

#include <stdio.h
#include <gsl/gsl_rng.h
#include <gsl/gsl_randist.h

int
main (void)
{
const gsl_rng_type * T;
gsl_rng * r;
int i, k;
int arrived = 0;
int accepted = 0;
int max_i = 10000000;
double u = 1.0;

gsl_rng_env_setup();
T = gsl_rng_default;
r = gsl_rng_alloc (T);

// In one iteration of this loop we simulate one hour.
for (i = 0; i < max_i; i++)
{
k = gsl_ran_poisson(r, u);

if (k)
{
arrived += k;
accepted += 1;
}
}

printf("u = %f\n", u);
printf("P(A|B) = %f\n", (double) accepted / (double) arrived);

return 0;
}

The function gsl_ran_poisson(r, u) returns a random integer according
to the Poisson distribution with the intensity u.


COMPARISON OF SIMULATION AND ANALYSIS
-------------------------------------

When you run the analysis, you'll get:

u = 1.000000
P(A|B) = 0.766988

while with the simulation:

u = 1.000000
P(A|B) = 0.632072


The two reported values of P(A|B) should be the same, and I'm unable
to explain why they are not. I would appreciate any advice.

Thanks for reading.


Best,
Irek


Back to top
ijs@iitis.gliwice.pl
science forum beginner


Joined: 29 Apr 2006
Posts: 10

PostPosted: Sun Apr 30, 2006 2:04 am    Post subject: Re: Poisson process: analysis and simulation disagree Reply with quote

I think I found what I do wrong in the analysis. When I compute
probabilities I disregard the fact that when there are k patients, the
probability of being accepted should contribute k times to the overall
average, while I was countig it only once by saying that:

P(A|B_k) = 1 / k.

Now my new analysis agrees with the simulation. I use the weighted
average. This time I don't go for conditional probabilities, but
instead say that I want to calculate an average probability that a
patient is accepted.

I derive the average probability from each of the following events:
- B_1: 1 patient arrived (probability f(1, u)),
- B_2: 2 patients arrived (probability f(2, u)),
..
..
..
- B_n: n patients arrived (probability f(n,u)).

I can't just average that way:

(1 + 1/2 + 1/3 + ... + 1/n) / n,

because events B_k are not equally probable. Moreover, in every event
there are different numbers of patients: 1, 2, ..., n. So I have the
following weight for event B_k.

w_k = f(k, u) * k

Finally, I calculate the weighted average this way:

sum of (w_k / k) for k = 1, ..., n
p = ----------------------------------
sum of (w_k) for k = 1, ..., n

My new program for the analysis:

#include <stdio.h>
#include <gsl/gsl_rng.h>
#include <gsl/gsl_randist.h>

int
main (void)
{
const gsl_rng_type * T;
gsl_rng * r;
int k, n = 10;
double u = 1.0;
double nominator = 0;
double denominator = 0;

gsl_rng_env_setup();
T = gsl_rng_default;
r = gsl_rng_alloc (T);

for (k = 1; k < n; k++)
{
// NOTE: k * (1.0 / (double)k) = 1
nominator += gsl_ran_poisson_pdf(k, u);
denominator += k * gsl_ran_poisson_pdf(k, u);
}

printf("u = %f\n", u);
printf("P = %f\n", nominator / denominator);
return 0;
}


Best,
Irek
Back to top
ijs@iitis.gliwice.pl
science forum beginner


Joined: 29 Apr 2006
Posts: 10

PostPosted: Sun Apr 30, 2006 2:16 am    Post subject: Re: Poisson process: analysis and simulation disagree Reply with quote

Dan and Mike,

Thank you for your posts. I fixed my analysis as I described in a post
a sent a while ago. I tried the simulation and analysis for various u
and they agree.

Thank you again!

Best,
Irek
Back to top
Ray Koopman
science forum Guru Wannabe


Joined: 25 Mar 2005
Posts: 216

PostPosted: Sun Apr 30, 2006 7:55 am    Post subject: Re: Poisson process: analysis and simulation disagree Reply with quote

Mike wrote:
Quote:
Can't follow your analysis and don't see why you define events A and B as
you do. You only need to focus on the hour in which YOU arrive in which
case the doctor will ALWAYS find at least one patient to choose from--YOU.
What is the point of defining Event A as "a patient is selected" and B as
"at least one patient arrives"?
P{A} = 1 and P{B} = 1.

Like Dan A, I agree with your simulation (not your analysis) and think of
the situation as follows:

Doctor opens his door once per hour and randomly selects one of the patients
waiting.
Patients not selected disappear instantly and new patients begin to arrive
and accumulate over the next hour until the doctor opens his door again to
make another selection. In this context, the Poisson distribution is
sensible for the number of patients (besides yourself) who have accumulated
in the waiting room over the hour, in that it allows the inter-arrival times
of patients to the waiting room to be exponentially distributed (average
rate of arrival = 1/u). (The exponential distribution is commonly used to
describe the time between successive arrivals of whatever and the Poisson
the number of such arrivals in a fixed period of time).

You, in effect, join a group whose size is (already) a Poisson variable.
Thus, if the "group" you join happens to be of size 0, your probability of
being selected is 100%. If the group you join is of size 1, your probability
of being selected is 50% and so on.

Your probability of being selected is exp(-u) times the summation of
(u^k)/(k+1)! from k=0 to infinity.
I didn't try to simply the series but I think it can be.
Anyway, a quick simulation produced the following (value of p agrees with
your simulated value when u = 1)
u= 0.01, p=0.9950;
u= 0.1, p=0.9516;
u= 0.5 p=0.7869;
u=1.0, p=0.6321;
u=2.0, p=0.4323;
u=3.0, p=0.3167;
u=10, p=0.1.
p seems to approach 1/u as u gets larger. Not sure why it doesn't approach
1/(u+1).

Let "A" be an arbitrarily specified person
in a population of n exchangeable people.

P(A sees the doctor | A is in the waiting room)

P(A is in the wr and sees the dr)
= ---------------------------------.
P(A is in the wr)


Denominator

= sum_{k=1}^n P(k people are in the wr) * P(A is one of those k)

= sum_{k=1}^n p(k) * ({n-1}C{k-1})/(nCk)

= sum_{k=1}^n p(k) * k/n.


Numerator

= sum_{k=1}^n P(k people are in the wr) * P(A is one of those k) / k

= sum_{k=1}^n p(k) / n.


Numerator
-----------
Denominator

sum_{k=1}^n p(k) / n
= ----------------------
sum_{k=1}^n p(k) * k/n

sum_{k=1}^n p(k) 1 - p(0)
= -------------------- = --------.
sum_{k=1}^n p(k) * k E(k)
Back to top
ijs@iitis.gliwice.pl
science forum beginner


Joined: 29 Apr 2006
Posts: 10

PostPosted: Mon May 01, 2006 7:08 pm    Post subject: Re: Poisson process: analysis and simulation disagree Reply with quote

Hi Mike,

Thanks for your posts.

Quote:
You only need to focus on the hour in which YOU arrive in which case
the doctor will ALWAYS find at least one patient to choose
from--YOU.

I should have made one thing clear in my previous post: I'm one of the
many patients, and I don't choose the hour at which I arrive. What I
wanted to say was that I'm interested in the probability that a
randomly chosen patient (just any waiting patient) from the waiting
room is accepted. So I'm interested in the average success
probability for patients in the waiting room. I didn't want to write
"what's the probability that any patient is selected," because it's 1:
a doctor will accept a patient if there is at least one waiting.
Therefore I wrote that I'm one of the patients in the waiting room. I
apologize for the confusion.


Best,
Irek
Back to top
Mike1170
science forum addict


Joined: 17 Sep 2005
Posts: 74

PostPosted: Tue May 02, 2006 1:54 am    Post subject: Re: Poisson process: analysis and simulation disagree Reply with quote

Irek wrote:
<ijs@iitis.gliwice.pl> wrote in message
news:1146510522.419423.326110@y43g2000cwc.googlegroups.com...
Quote:
I should have made one thing clear in my previous post: I'm one of the
many patients, and I don't choose the hour at which I arrive. What I
wanted to say was that I'm interested in the probability that a
randomly chosen patient (just any waiting patient) from the waiting
room is accepted. So I'm interested in the average success
probability for patients in the waiting room. I didn't want to write
"what's the probability that any patient is selected," because it's 1:
a doctor will accept a patient if there is at least one waiting.
Therefore I wrote that I'm one of the patients in the waiting room. I
apologize for the confusion.

In that case, I think the approach suggested by Ray Koopman is the one to
use but what about the inconsistency between a distribution (Poisson)
defined on an infinite domain but a finite population from which to create
it. Nevermind. I suppose for large n and small u, the problem is only
theoretical. But then again, the problem seems to be equivalent to the
following (which makes the whole problem trivial):

A. If k=0 (no patients in waiting room), then the probability of the
doctor selecting any given person = 0. This situation occurs with
probability exp(-u). P{k=0} for Poisson
B. If k>0, then probability of doctor selecting any given person = 1/n.
See below. This situation occurs with probability 1-exp(-u).

Using the approach suggested by Ray K (a hypothetical population of n
patients from which to choose the k waiting-room patients), the probability
that the doctor selects a given patient from the k patients in the waiting
room (which is itself a random sample of k patients chosen from n patients)
is equivalent to the probability that the doctor selects a given patient
directly from the entire population of n patients.
The prob given patient is in waiting room = ( (n-1)C(k-1) x 1C1 ) / nCk =
k/n
The prob given patient is selected by Dr. = ( (k-1)C0 ) x 1C1 ) / kC1 =
1/k
k / n x 1 / k = 1 / n, right?
Back to top
ijs@iitis.gliwice.pl
science forum beginner


Joined: 29 Apr 2006
Posts: 10

PostPosted: Wed May 03, 2006 1:25 am    Post subject: Re: Poisson process: analysis and simulation disagree Reply with quote

Ray,

Thank you for your post. Your analysis is interesting, especially how
you got the solution for the infinite population based on the analysis
for a finite population of n elements. That was enlightening! Thanks
again.

Best,
Irek
Back to top
danheyman@yahoo.com
science forum beginner


Joined: 18 Jul 2005
Posts: 33

PostPosted: Sat May 06, 2006 12:51 am    Post subject: Re: Poisson process: analysis and simulation disagree Reply with quote

This analysis is correct, and the weighting is the heart of the issue.
The fact the the probability that a patient is in a batch of exactly k
patients is not the probability that a batch has k patients (let's
denote that that by p(k)) but is k*p(k)/u, where u is the mean batch
size, is known as the "inspection paradox". This holds for all
distributions on the non-negative integers, so the answer to your
problem is [1-p(0)]/u which is 1-exp(-1) for the Poisson dst. with mean
1.

Dan Heyman
Back to top
Google

Back to top
Display posts from previous:   
Post new topic   Reply to topic Page 1 of 1 [15 Posts] View previous topic :: View next topic
The time now is Fri Oct 20, 2017 7:45 pm | All times are GMT
Forum index » Science and Technology » Math » Probability
Jump to:  

Similar Topics
Topic Author Forum Replies Last Post
No new posts electroplating organic chemicals analysis a.1998.a@gmail.com Chem 0 Fri Jul 21, 2006 11:49 am
No new posts reference books for formulaes in stochastic analysis and ... Michael11 Math 0 Thu Jul 20, 2006 12:38 am
No new posts Functional analysis? e_hobsbawm@hotmail.com Math 1 Wed Jul 19, 2006 1:55 pm
No new posts Is there a way to write out the process of the cumulative... Michael11 Math 1 Wed Jul 19, 2006 7:16 am
No new posts Trying to Understand the Set Theory of Non-Standard Analysis Math1723 Math 1 Wed Jul 19, 2006 3:52 am

Copyright © 2004-2005 DeniX Solutions SRL
Other DeniX Solutions sites: Electronics forum |  Medicine forum |  Unix/Linux blog |  Unix/Linux documentation |  Unix/Linux forums  |  send newsletters
 


Powered by phpBB © 2001, 2005 phpBB Group
[ Time: 0.0592s ][ Queries: 16 (0.0078s) ][ GZIP on - Debug on ]