Poisson distribution: calculate lambda

[drfg2008]

A question that is perhaps not adequate for this platform.

I try to  dissolve the Poisson distribution for Lambda. Unfortunately I do not not have access to Mathematica right now. Does anyone know the formula by chance?

The problem:

I have  p, k.  What is lambda?

Even better would be a solution to find lambda in the following case:

I have:
the sum of the probabilities p (i) up to a certain k
and search Lambda .
(example: I have as the sum of p(i) = 50% and a k = 5. What would be Lambda ? )


As I said, the question is probably not really adequate here. But maybe someone knows a solution. Maybe even with SPSS.

Frank

[Art Kendall]

an approximation can be found with the syntax below.
The first file steps by integers, using the results of that you would guess the ballpark looping to find more significant digits.

There are several ways this can be refined depending how you are going to use it.
E.g., very fine steps, but only selecting for printing or matching to another file the lambda with the smallest abs(fudge), etc.


Art Kendall
Social Research Consultants

new file.
input program.
   loop lambda = 1 to 25.
      compute k =5.
      compute mycumprob = .50.
      compute cumprob = cdf.poisson(k,lambda).
      end case.
   end loop.
   end file.
end input program.
compute fudge = mycumprob-cumprob.
compute constant = 1.
formats lambda (f5.3) k constant (f2) mycumprob cumprob fudge (f5.4).
dataset name run1.
aggregate outfile= * mode=addvariables /break = constant
   /minfudge = min(fudge).
list variables = lambda fudge.
new file.
input program.
   loop lambda = 5.3 to 5.9 by .025.
      compute k =5.
      compute mycumprob = .50.
      compute cumprob = cdf.poisson(k,lambda).
      end case.
   end loop.
   end file.
end input program.
compute fudge = mycumprob-cumprob.
compute constant = 1.
formats lambda (f5.3) k constant (f2) mycumprob cumprob fudge (f5.4).
aggregate outfile= * mode=addvariables /break = constant
   /minfudge = min(fudge).
dataset name run2.
list variables = lambda fudge.


On 11/22/2011 5:36 AM, drfg2008 wrote:

A question that is perhaps not adequate for this platform.

I try to  dissolve the Poisson distribution for Lambda. Unfortunately I do
not not have access to Mathematica right now. Does anyone know the formula
by chance?

The problem:

I have  p, k.  What is lambda?

Even better would be a solution to find lambda in the following case:

I have:
the sum of the probabilities p (i) up to a certain k
and search Lambda .
(example: I have as the sum of p(i) = 50% and a k = 5. What would be Lambda
? )


As I said, the question is probably not really adequate here. But maybe
someone knows a solution. Maybe even with SPSS.

Frank

-----
Dr. Frank Gaeth
FU-Berlin

--
View this message in context: http://spssx-discussion.1045642.n5.nabble.com/Poisson-distribution-calculate-lambda-tp5013192p5013192.html
Sent from the SPSSX Discussion mailing list archive at Nabble.com.

=====================
To manage your subscription to SPSSX-L, send a message to
[hidden email] (not to SPSSX-L), with no body text except the
command. To leave the list, send the command
SIGNOFF SPSSX-L
For a list of commands to manage subscriptions, send the command
INFO REFCARD

===================== To manage your subscription to SPSSX-L, send a message to [hidden email] (not to SPSSX-L), with no body text except the command. To leave the list, send the command SIGNOFF SPSSX-L For a list of commands to manage subscriptions, send the command INFO REFCARD

[David Marso]

Somewhat of a refinement of Art's post:
HTH, David
---
DEFINE POISP (K !TOKENS(1) / P !TOKENS(1) / TOL !TOKENS(1)).

NEW FILE.
INPUT PROGRAM.
COMPUTE k =!K .
COMPUTE P = !P .
LOOP lambda = 1 TO 25.
+    COMPUTE cumprob = CDF.POISON(k,lambda).
+    LEAVE k P.
+    END CASE.
END LOOP.
END FILE.
END INPUT PROGRAM.

COMPUTE fudge = ABS(P-cumprob).

COMPUTE delta=(lambda-lag(lambda))/20.
SORT CASES BY FUDGE.
SELECT IF $CASENUM=1.
LOOP.
+  COMPUTE MIN=1.
+  LOOP Lambda2=lambda-10*delta TO lambda+10*delta BY delta .
+    COMPUTE cumprob = CDF.POISON(k,Lambda2).
+    COMPUTE fudge=ABS(P-cumprob).
+    DO IF fudge < MIN.
+      COMPUTE L=Lambda2.
+      COMPUTE MIN=fudge.
+    END IF.
+  END LOOP.
+  COMPUTE LAMBDA=L.
+  COMPUTE delta=delta/20.
END LOOP IF MIN < !TOL .
PRINT / "K=" K /"P=" P / "LAMBDA=" LAMBDA (F12.8) /"ERROR=" MIN (F12.8).
EXE.
!ENDDEFINE.
POISP K=5 P =.5 TOL=.00000001 .
POISP K=6 P =.5 TOL=.00000001 .
POISP K=7 P =.5 TOL=.00000001 .
POISP K=8 P =.5 TOL=.00000001 .

Art Kendall wrote

an approximation can be found with the syntax below.
      The first file steps by integers, using the results of that you
      would guess the ballpark looping to find more significant digits.
     
      There are several ways this can be refined depending how you are
      going to use it.
      E.g., very fine steps, but only selecting for printing or matching
      to another file the lambda with the smallest abs(fudge), etc.
     
     
      Art Kendall
      Social Research Consultants
     
      new file.
      input program.
         loop lambda = 1 to 25.
            compute k =5.
            compute mycumprob = .50.
            compute cumprob = cdf.poisson(k,lambda).
            end case.
         end loop.
         end file.
      end input program.
      compute fudge = mycumprob-cumprob.
      compute constant = 1.
      formats lambda (f5.3) k constant (f2) mycumprob cumprob fudge
      (f5.4).
      dataset name run1.
      aggregate outfile= * mode=addvariables /break = constant
         /minfudge = min(fudge).
      list variables = lambda fudge.
      new file.
      input program.
         loop lambda = 5.3 to 5.9 by .025.
            compute k =5.
            compute mycumprob = .50.
            compute cumprob = cdf.poisson(k,lambda).
            end case.
         end loop.
         end file.
      end input program.
      compute fudge = mycumprob-cumprob.
      compute constant = 1.
      formats lambda (f5.3) k constant (f2) mycumprob cumprob fudge
      (f5.4).
      aggregate outfile= * mode=addvariables /break = constant
         /minfudge = min(fudge).
      dataset name run2.
      list variables = lambda fudge.
   
   
    On 11/22/2011 5:36 AM, drfg2008 wrote:
   
      A question that is perhaps not adequate for this platform.

I try to  dissolve the Poisson distribution for Lambda. Unfortunately I do
not not have access to Mathematica right now. Does anyone know the formula
by chance?

The problem:

I have  p, k.  What is lambda?

Even better would be a solution to find lambda in the following case:

I have:
the sum of the probabilities p (i) up to a certain k
and search Lambda .
(example: I have as the sum of p(i) = 50% and a k = 5. What would be Lambda
? )


As I said, the question is probably not really adequate here. But maybe
someone knows a solution. Maybe even with SPSS.

Frank

-----
Dr. Frank Gaeth
FU-Berlin

--
View this message in context: http://spssx-discussion.1045642.n5.nabble.com/Poisson-distribution-calculate-lambda-tp5013192p5013192.html 
Sent from the SPSSX Discussion mailing list archive at Nabble.com.

=====================
To manage your subscription to SPSSX-L, send a message to
[hidden email]  (not to SPSSX-L), with no body text except the
command. To leave the list, send the command
SIGNOFF SPSSX-L
For a list of commands to manage subscriptions, send the command
INFO REFCARD


   
 


=====================
To manage your subscription to SPSSX-L, send a message to
[hidden email] (not to SPSSX-L), with no body text except the
command. To leave the list, send the command
SIGNOFF SPSSX-L
For a list of commands to manage subscriptions, send the command
INFO REFCARD

[drfg2008]

Thank you, as always a great help.

My problem is that I have tens of thousands of cases for which I have to calculate lambda. To calculate each individual case by iterations will mean a great effort. So my question is whether there is a possibility with SPSS to look into external tables (matrix) to a certain criterion (here: the probability p) and to pick up a value in the same row in order to read out the lambda value.

Excel has a function like this:

 VLOOKUP (E16, E9: F12, 2, TRUE)


Frank

[Art Kendall]

perhaps if you describe the overall problem list members could suggest ways to try to answer the underlying questions.
Are the things that are common to the problem etc.?

Art Kendall
 

On 11/24/2011 5:18 AM, drfg2008 wrote:

Thank you, as always a great help.

My problem is that I have tens of thousands of cases for which I have to
calculate lambda. To calculate each individual case by iterations will mean
a great effort. So my question is whether there is a possibility with SPSS
to look into external tables (matrix) to a certain criterion (here: the
probability p) and to pick up a value in the same row in order to read out
the lambda value.

Excel has a function like this:

 VLOOKUP (E16, E9: F12, 2, TRUE)


Frank

-----
Dr. Frank Gaeth
FU-Berlin

--
View this message in context: http://spssx-discussion.1045642.n5.nabble.com/Poisson-distribution-calculate-lambda-tp5013192p5019778.html
Sent from the SPSSX Discussion mailing list archive at Nabble.com.

=====================
To manage your subscription to SPSSX-L, send a message to
[hidden email] (not to SPSSX-L), with no body text except the
command. To leave the list, send the command
SIGNOFF SPSSX-L
For a list of commands to manage subscriptions, send the command
INFO REFCARD

===================== To manage your subscription to SPSSX-L, send a message to [hidden email] (not to SPSSX-L), with no body text except the command. To leave the list, send the command SIGNOFF SPSSX-L For a list of commands to manage subscriptions, send the command INFO REFCARD

[drfg2008]

Art,

I didn't quite understand your mail.

As I understand the syntax (above) this solves the problem to compute lambda in one case (maybe I'm wrong?).

Now, since I have a table of more than ten thousand cases, where I have to compute a lambda for each case, I thought of creating a 'lambda-table' with three columns: p / k  / lambda and 'cross'-check it with my ten thousand cases. With excel this procedure would be possible with the command VLOOKUP (E16, E9: F12, 2, TRUE) ( German: sverweis).

My question is if in SPSS there is a similar function that would look up the value out of a list.

So, for example, in case no. 123 the p-value might be 0,5 an dk = 5. In excel one could 'match' the two tables via VLOOKUP.

Probably this is a too excel-like thinking. Sorry for that.

Frank

[Jon K Peck]

There is a function named - surprise - vlookup in the extendedTransforms.py module available on the SPSS Community site.

Here is a usage example from the doc in the module
# Example using the SPSSINC TRANS extension command
#* The lookup table.
#data list free/ value(F8.0) akey(A1).
#begin data
#10 'a'
#20 'b'
#100 'z'
#end data.
#dataset name lookup.

#* The main dataset.
#data list free/x(f8.0) y(A2).
#begin data
#1 'a'
#2 'b'
#5 'a '
#10 ''
#1 'b'
#end data.
#dataset name main.
#dataset activate main.

#spssinc trans result = resultcodealpha
#/initial "extendedTransforms.vlookup('akey', 'value', 'lookup')"
#/formula func(y).

HTH,

Jon Peck (no "h") aka Kim
Senior Software Engineer, IBM
[hidden email]
new phone: [phone]




From:        drfg2008 <[hidden email]>
To:        [hidden email]
Date:        11/24/2011 03:21 AM
Subject:        Re: [SPSSX-L] Poisson distribution: calculate lambda
Sent by:        "SPSSX(r) Discussion" <[hidden email]>

---

Thank you, as always a great help.

My problem is that I have tens of thousands of cases for which I have to
calculate lambda. To calculate each individual case by iterations will mean
a great effort. So my question is whether there is a possibility with SPSS
to look into external tables (matrix) to a certain criterion (here: the
probability p) and to pick up a value in the same row in order to read out
the lambda value.

Excel has a function like this:

VLOOKUP (E16, E9: F12, 2, TRUE)


Frank

-----
Dr. Frank Gaeth
FU-Berlin

--
View this message in context: http://spssx-discussion.1045642.n5.nabble.com/Poisson-distribution-calculate-lambda-tp5013192p5019778.html
Sent from the SPSSX Discussion mailing list archive at Nabble.com.

=====================
To manage your subscription to SPSSX-L, send a message to
[hidden email] (not to SPSSX-L), with no body text except the
command. To leave the list, send the command
SIGNOFF SPSSX-L
For a list of commands to manage subscriptions, send the command
INFO REFCARD

[drfg2008]

Great. I knew something like that must exist.

Frank

PS Sorry if my descriptions were a little imprecise and too excel-like.

[David Marso]

"To calculate each individual case by iterations will mean a great effort. "
Would probably take only a few seconds.  May need to reorganize my code a little but should be trivial.
How are you intending to generate the lookup table in the first place?
--

drfg2008 wrote

Thank you, as always a great help.

My problem is that I have tens of thousands of cases for which I have to calculate lambda. To calculate each individual case by iterations will mean a great effort. So my question is whether there is a possibility with SPSS to look into external tables (matrix) to a certain criterion (here: the probability p) and to pick up a value in the same row in order to read out the lambda value.

Excel has a function like this:

 VLOOKUP (E16, E9: F12, 2, TRUE)


Frank

[drfg2008]

but should be trivial ...

well, for you. I'm sure.

This is how the data looks like:

input program.
loop gameNo =1 to 10000 by 1.
end case.
end loop.
end file.
end input program.
EXECUTE.

COMPUTE Team1=RV.UNIFORM(1,6)/100.
COMPUTE k = 5.
COMPUTE lambda = 123.
EXECUTE.


Team1 column consists of probabilities. The p is cumulative p.

For example, if Team1 had p=0.548510246  and k = 5, lambda would be 5.3863.

[Art Kendall]

Please tell us what you are trying to simulate.


Should the compute for team1 find integers?
If it is that would simplify calculating a lookup table.

In a dice simulation it would be
COMPUTE die1=RND(RV.UNIFORM(.5,6.5)).

Art Kendall
Social Research Consultants


On 11/25/2011 5:17 AM, drfg2008 wrote:

> /but should be trivial .../
>
> well, for you. I'm sure.
>
> This is how the data looks like:
>
> input program.
> loop gameNo =1 to 10000 by 1.
> end case.
> end loop.
> end file.
> end input program.
> EXECUTE.
>
> COMPUTE Team1=RV.UNIFORM(1,6)/100.
> COMPUTE k = 5.
> COMPUTE lambda = 123.
> EXECUTE.
>
>
> Team1 column consists of probabilities. The p is cumulative p.
>
> For example, if Team1 had p=0.548510246  and k = 5, lambda would be 5.3863.
>
>
>
> -----
> Dr. Frank Gaeth
> FU-Berlin
>
> --
> View this message in context: http://spssx-discussion.1045642.n5.nabble.com/Poisson-distribution-calculate-lambda-tp5013192p5022448.html
> Sent from the SPSSX Discussion mailing list archive at Nabble.com.
>
> =====================
> To manage your subscription to SPSSX-L, send a message to
> [hidden email] (not to SPSSX-L), with no body text except the
> command. To leave the list, send the command
> SIGNOFF SPSSX-L
> For a list of commands to manage subscriptions, send the command
> INFO REFCARD
>

=====================
To manage your subscription to SPSSX-L, send a message to
[hidden email] (not to SPSSX-L), with no body text except the
command. To leave the list, send the command
SIGNOFF SPSSX-L
For a list of commands to manage subscriptions, send the command
INFO REFCARD

[David Marso]

What is the range of values for K?
What have you tried so far to apply what I posted recently to your N=10000?
What impedes your approach to a solution.
I think a VLOOKUP idea generally sucks because where do you get the table?  
Probably too granular (inprecise) and not particularly elegant by any stretch of the imagination!!!

Hint:  Rather than INPUT PROGRAM for generation of cases in the first part of my solution, build it horizontally, locate initial integer value then apply the iterative code to that value.  Should be utterly trivial since all the hard work has already been done.  You simply need to study it and apply it.

drfg2008 wrote

but should be trivial ...

well, for you. I'm sure.

This is how the data looks like:

input program.
loop gameNo =1 to 10000 by 1.
end case.
end loop.
end file.
end input program.
EXECUTE.

COMPUTE Team1=RV.UNIFORM(1,6)/100.
COMPUTE k = 5.
COMPUTE lambda = 123.
EXECUTE.


Team1 column consists of probabilities. The p is cumulative p.

For example, if Team1 had p=0.548510246  and k = 5, lambda would be 5.3863.

[David Marso]

Try this on for size.  Like I said,  utterly trivial ;-)
--
*Data Simulation *.
input program.
loop gameNo =1 to 10000.
COMPUTE T_Prob=UNIFORM(1).
COMPUTE k = TRUNC(UNIFORM(10))+3.
end case.
end loop.
end file.
end input program.
* Check for Frank's example *.
DO IF $CASENUM=1.
COMPUTE K=5.
COMPUTE T_Prob=0.548510246.
END IF.
EXE.


*Estimate Lambda -> Lambda_F *.
NUMERIC LAMBDA_F ABS_ERR.
VECTOR #P(30).
LOOP #LAMBDA0=1 TO 30.
COMPUTE #P(#LAMBDA0)=CDF.POISON(k,#lambda0).
END LOOP If #P(#LAMBDA0) < T_Prob.
COMPUTE #delta=1/10.
COMPUTE #LAMBDA=#LAMBDA0-1.

***LOOP /* This seems to hang sometimes ;-((( */ So use LOOP To 5  instead.
LOOP #=1 To 5.
+  COMPUTE #MIN=1.
+  LOOP #Lambda2=#lambda-10*#delta TO #lambda+10*#delta BY #delta .
+    COMPUTE #fudge=ABS(T_PROB-CDF.POISON(k,#Lambda2)).
+    DO IF #fudge < #MIN.
+      COMPUTE #L=#Lambda2.
+      COMPUTE #MIN=#fudge.
+    END IF.
+  END LOOP.
+  COMPUTE #LAMBDA=#L.
+  COMPUTE #delta=#delta/20.
END LOOP.
***END LOOP IF #MIN< .00001  /* This seems to hang sometimes ;-((( */.
COMPUTE LAMBDA_F=#LAMBDA.
COMPUTE ABS_ERR=#MIN.
EXE.
* SORT CASES BY ABS_Err(D).

David Marso wrote

What is the range of values for K?
What have you tried so far to apply what I posted recently to your N=10000?
What impedes your approach to a solution.
I think a VLOOKUP idea generally sucks because where do you get the table?  
Probably too granular (inprecise) and not particularly elegant by any stretch of the imagination!!!

Hint:  Rather than INPUT PROGRAM for generation of cases in the first part of my solution, build it horizontally, locate initial integer value then apply the iterative code to that value.  Should be utterly trivial since all the hard work has already been done.  You simply need to study it and apply it.

drfg2008 wrote

but should be trivial ...

well, for you. I'm sure.

This is how the data looks like:

input program.
loop gameNo =1 to 10000 by 1.
end case.
end loop.
end file.
end input program.
EXECUTE.

COMPUTE Team1=RV.UNIFORM(1,6)/100.
COMPUTE k = 5.
COMPUTE lambda = 123.
EXECUTE.


Team1 column consists of probabilities. The p is cumulative p.

For example, if Team1 had p=0.548510246  and k = 5, lambda would be 5.3863.

[David Marso]

SCRAP that last post ;-(  I screwed the pooch on the upper bound for the inner iteration.  Wondering why it didn't converge *DOH... Channeling Homer Simpson*.  
Try out this lean, mean revision ;-)
9 lines of code (Using T_Prob for target probability, K for K)-

*Estimate Lambda -> Lambda_F *.
COMPUTE #delta=1.
LOOP LAMBDA_F=1 TO 50 BY #delta.
END LOOP IF CDF.POISSON(k,LAMBDA_F)<= T_Prob.
LOOP .
+  COMPUTE #delta=#delta/10.
+  LOOP Lambda_F=MAX(.000001,LAMBDA_F-10*#delta) TO LAMBDA_F+10*#delta BY #delta .
+  COMPUTE ERROR=CDF.POISSON(k,Lambda_F) - T_Prob.
+  END LOOP IF ERROR LE 0.
END LOOP IF ABS(ERROR) <=.00000001 .
EXE.

[David Marso]

This version is more general and avoids a few problems (neg Lambda and restrictions on large K).
VERY ACCURATE!!!!

*Estimate Lambda -> Lambda_F *.
SET MXLOOPS=1000.
COMPUTE #delta=1.
COMPUTE LAMBDA_F=0.
LOOP .
+  LOOP .
+    COMPUTE LAMBDA_F=MAX(.00001,LAMBDA_F+#delta).
+  END LOOP IF (CDF.POISSON(k,LAMBDA_F) - T_Prob) LE 0.
+  COMPUTE LAMBDA_F=MAX(.00001,LAMBDA_F-#delta).
+  COMPUTE #delta=#delta/10.
END LOOP IF ABS(CDF.POISSON(k,LAMBDA_F) - T_Prob  ) <=.0000000001 .
COMPUTE ABS_ERR=ABS(CDF.POISSON(k,LAMBDA_F) - T_Prob  ) .
EXE.  

-----

David Marso wrote

SCRAP that last post ;-(  I screwed the pooch on the upper bound for the inner iteration.  Wondering why it didn't converge *DOH... Channeling Homer Simpson*.  
Try out this lean, mean revision ;-)
9 lines of code (Using T_Prob for target probability, K for K)-

*Estimate Lambda -> Lambda_F *.
COMPUTE #delta=1.
LOOP LAMBDA_F=1 TO 50 BY #delta.
END LOOP IF CDF.POISSON(k,LAMBDA_F)<= T_Prob.
LOOP .
+  COMPUTE #delta=#delta/10.
+  LOOP Lambda_F=MAX(.000001,LAMBDA_F-10*#delta) TO LAMBDA_F+10*#delta BY #delta .
+  COMPUTE ERROR=CDF.POISSON(k,Lambda_F) - T_Prob.
+  END LOOP IF ERROR LE 0.
END LOOP IF ABS(ERROR) <=.00000001 .
EXE.

[David Marso]

Drastically optimized version!!!!!!!!!  I tested with N=1000000.
10000 runs in about 5 seconds on my ancient MacBookPro (circa 2007).
Windoze Vista.
------
*Estimate Lambda *.
SET MXLOOPS=100.
COMPUTE #delta=1.
LOOP .
+  LOOP .
+    COMPUTE LAMBDA = MAX(.00001,LAMBDA + #delta).
+  END LOOP IF CDF.POISSON(k,LAMBDA) LE T_Prob.
+  COMPUTE LAMBDA = LAMBDA-#delta.
+  COMPUTE #delta=#delta/3 .
END LOOP IF ABS(CDF.POISSON(k,LAMBDA) - T_Prob  ) LE .0000000001 .
EXE.

[drfg2008]

That's it. But far from easy. At least for me.
Thanks.

Frank