2
$\begingroup$

How to plot the following in Mathematica:

$$P = \sum_{k=C+1}^{\infty}\frac{e^{-C^\gamma} (C^\gamma)^k}{k!}$$

I want to plot $P$ versus $C \log{C}$ for a given $\gamma$ with logarithmic scale on Y-axis. In MATLAB, I'd calculate it as follows:

gamma = 0.75;
C = 1:1:20;
lambda = C .^ gamma;
p = 1 - poisscdf(C,lambda);
semilogy(C .* log(C), p);
$\endgroup$
  • $\begingroup$ You probably could have worked this out by searching in the documentation for "poisson distribution", yielding PoissonDistribution, and "log plot", yielding a range of hits including ListLogPlot, LogPlot etc. $\endgroup$ – Verbeia May 3 '12 at 4:38
11
$\begingroup$

It is worth noting two characteristic features of the question:

  1. The sum is a complementary cumulative distribution function for a Poisson distribution: it's built in to Mathematica and needn't be computed explicitly.

  2. The use of $C+1$ as a starting index in the sum, as well as the expression 1:1:20 in the code, indicate $C$ is considered an integer: this needs to be a discrete plot.

It can also help to draw clear parallels between the ΜATLAB approach and an idiomatic Mathematica approach. How about this?

Module[{c = Range[20], γ = 0.75, p},
 p = {# Log[#], 1 - CDF[PoissonDistribution[#^γ], #]} & /@ c;
 ListLogPlot[p]
]

Log linear plot

| improve this answer | |
$\endgroup$
  • $\begingroup$ Much easier now. But, what about the logarithmic scale on the y-axis? I mean the ticks! $\endgroup$ – Osama Gamal May 2 '12 at 14:19
  • 4
    $\begingroup$ @Osama They are logarithmic! $\endgroup$ – whuber May 2 '12 at 19:14
6
$\begingroup$

Here is what I would do:

p[c_, x_] = (Sum[Exp[-x] x^k/k!, {k, c + 1, Infinity}] // 
   FullSimplify)

With[{gamma = 0.75}, 
 ParametricPlot[{c Log[c], p[c, c^gamma]}, {c, 1, 20}, 
  AspectRatio -> 1]]

enter image description here

| improve this answer | |
$\endgroup$
  • $\begingroup$ How to make it logarithmic scale on y axes as in the Matlab code? I used Log[p[c, c^gamma]] but still can't edit the axes thing!! $\endgroup$ – Osama Gamal May 2 '12 at 13:32
  • $\begingroup$ cValues = Range[1, 20, 0.05]; ListLogPlot[{# Log[#], p[#, #^0.75]} & /@ cValues, PlotRange -> All] ? $\endgroup$ – b.gates.you.know.what May 2 '12 at 13:41
  • $\begingroup$ Why things is much complicated in Mathematica than in Matlab! $\endgroup$ – Osama Gamal May 2 '12 at 14:12
  • 1
    $\begingroup$ Here's a more compact definition of p[c, x]: p[c_, x_] := GammaRegularized[c + 1, 0, x]. I don't understand why FullSimplify[] doesn't use GammaRegularized[] more frequently myself... $\endgroup$ – J. M.'s technical difficulties May 2 '12 at 14:25
3
$\begingroup$

This is only a minor variant on the existing answers, but it will hopefully be clearer to a novice user:

pp[c_?Positive, γ_] := {c Log[c], 1 - CDF[PoissonDistribution[c^γ], c]}
ListLogPlot[Table[pp[x, 0.75], {x, 1, 40}]]

enter image description here

But notice:

ListLogPlot[Table[pp[x, 0.75], {x, 1, 40, 0.02}], Joined -> True]

enter image description here

The function is defined for non-integer values of $c$, but includes a Floor expression.

Assuming[c > 0, Simplify[1 - CDF[PoissonDistribution[c^γ], c]]]
(*  ==> 1 - GammaRegularized[1 + Floor[c], c^γ]  *)

As J.M. pointed out in comments, the expression that results from this simplification can be further simplified to GammaRegularized[1 + Floor[c], 0, c^γ], but the presence of Floor[] in the result is clear even in the less-simplified version.

| improve this answer | |
$\endgroup$
  • 1
    $\begingroup$ Again: 1 - GammaRegularized[1 + Floor[c], c^γ] can be simplified further to GammaRegularized[1 + Floor[c], 0, c^γ]. Numerically, it's a good idea to be using the three-argument form GammaRegularized[] for situations like this... $\endgroup$ – J. M.'s technical difficulties May 3 '12 at 6:43

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.