Skip to content

Instantly share code, notes, and snippets.

@VincentTam
Last active August 29, 2015 14:18
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save VincentTam/c27c38c49fe8de17c815 to your computer and use it in GitHub Desktop.
Save VincentTam/c27c38c49fe8de17c815 to your computer and use it in GitHub Desktop.
Calculate P(A > B) and P(A = B), where A and B are independent Poisson variables.
#include <iostream>
#include <cmath>
using namespace std;
double pXy(double x, double y, int N);
double pxy(double x, double y, int N);
int main(void) {
double pAb,paB,pab,a,b;
int N;
cout << "Assume that Poisson r.v. A and B are indepedent" << endl;
cout << "Parameter for A: ";
cin >> a;
cout << "Parameter for B: ";
cin >> b;
cout << "Number of terms to be added (100 <= N <= 1000): ";
cin >> N;
pAb = pXy(a,b,N);
paB = pXy(b,a,N);
pab = pxy(a,b,N);
cout << "P(A > B) = " << pAb << ", P(A < B) = " << paB <<
", P(A = B) = " << pab << endl;
}
/* P(X > Y) */
double pXy(double x, double y, int N) {
double ans = 0;
for (int i = 1; i <= N; i++) {
for (int j = 0; j < i ; j++) {
double term = 1;
for (int k = 1; k <= i; k++)
term *= x / k;
for (int k = 1; k <= j; k++)
term *= y / k;
ans += term;
}
}
return ans * exp(-x - y);
}
/* P(X = Y) */
double pxy(double x, double y, int N) {
double ans = 0;
for (int i = 0; i <= N; i++) {
double term = 1;
for (int k = 1; k <= i; k++)
term *= x / k * y / k;
ans += term;
}
return ans * exp(-x - y);
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment