You should not worry about overdispersion in such a model. Poisson regression can be used because the likelihood function (as a function of parameters) is the same as what you would have if data was poisson distributed. But the original data is time-to-event data, and only because of a kind of mathematical luck, the parameters can be estimated by use of Poisson regression. That means, you can make statistical inference about parameters with use of Poisson-regression, but you can not use the statistics about how good your data is Poisson distributed.   For example, if you have N observations with time to event and you will estimate the intensity, then you put N on the left side (which is of course not Poisson distributed) and the log of the total observation time in offset. It would of course then be meaningless to test how good N is Poisson distributed since it is a fixed number. That example is ofcourse rather extreme, but it illustrate the principle. ... View more

Hello,     Here is some methods you can use to validate the assumption in the PH-model, but they are all cumbersome when there is many independent variables. You can divide the time axis into pieces (two pieces for instance). Then check interaction between the time and the variable. Its a good idea to aggregate your data on riskset and then use the type3-option to check interaction (see this link). Otherwise you should define the timedependence inside phreg  (which increase calculation alot). You will get a p-value for each of you independent variables.   An other way to check the assumption is by using the assess statement. Here you check that the score-process behaves as it should if the assumption is fulfilled. You will get a plot of the score-process for each parameter in the model. The scoreprocess should look like the simulatied curves (see example below). Calculation time is proportional to the square of events, so if you have a large datasets this will no be an option.   data simulation; do i=1 to 500; t=rand('exponential',10); group=rand('bernoulli',0.5); output; end; run; ods html; ods graphics; proc phreg data=simulation; class group/param=glm; model t=group; assess ph/npaths=100 resample=1; run; ods html close;   It is worthwhile to mention the that you can plot log(-log(S1(t)), where S(t) is kaplan meier plots. The plot you get should be constant if the assumption is fulfilled. This method only work with one variable at a time. Here is an example:   　 data simulation; do i=1 to 20000; group=rand('bernoulli',0.5); t=rand('exponential',10*exp(-0.5*group)); output; end; run; proc sort data=simulation; by group; run; proc phreg data=simulation; class group/param=glm; model t=; baseline out=baseline survival=survival/ method=pl; by group; run; data dif; set baseline(in=a0 where=(group=0) rename=(survival=s0)) baseline(in=a1 where=(group=1) rename=(survival=s1)); by t; retain lasts0 lasts1; if a0 then s1=lasts1; else if a1 then s0=lasts0; dif=log(-log(s0))-log(-log(s1)); if a0 then lasts0=s0; else if a1 then lasts1=s1; run; 　 symbol i=join v=none; goptions device=win; proc gplot data=dif; plot dif*t; run;   I will mention that I find it funny that there are so much focus on whether there is proportional hazards. Why is this interaction between some indendent variables and time more important than interaction between the independent variabels. That is because Cox's proporportional hazard model is named after this assumption! For instance, the same assumption is used in Poisson regression model, but reviewers rarely ask to check the proportional hazard if Poisson regression is used. ... View more

Or try this data elipse;   pi=4*atan(1);   do i=0 to 256 ;     x=cos(i*2*pi/256);     y=2*sin(i*2*pi/256); output;   end; run; data _NULL_;   set elipse end=eof;   retain area lastx lasty;   if _N_=1 then area=0;   else do;     area+((lastx-x)*(y+lasty)/2);   end;   lastx=x;lasty=y;   if eof then put area=; run; ... View more

You get the confidence intervals of the Log(IRR) in the parmest dataset. They can just be expontiated and you will have confidence intervals of IRR. If you want the confidence intervals calculated within proc genmod you can use the estimate statement. data abc;   do i=1 to 100;     group=int(5*ranuni(-1));     y=rand('poisson',exp(group));     output;   end; run; ods output parameterestimates=parmest; proc genmod data=abc;   class group;   model y=group/dist=poisson link=log; run; data parmest;   set parmest;   text=put(exp(estimate),5.2)||'('||put(exp(lowerwaldcl),5.2)||','||put(exp(upperwaldcl),5.2)||')';   put parameter= @25 level1= text=; run; *and if you want the intervals calculated by genmod:; ods output estimates=estimates; proc genmod data=abc;   class group;   model y=group/dist=poisson link=log;   estimate 'group 0 vs 4' group 1 0 0 0 -1/exp;   estimate 'group 1 vs 4' group 0 1 0 0 -1/exp;   estimate 'group 2 vs 4' group 0 0 1 0 -1/exp;   estimate 'group 3 vs 4' group 0 0 0 1 -1/exp; run; ... View more

It is likely that the difference is due to the covariance structure. You used a unstructured covariance matrix, but I would guess that STATA use the varicance-component structure as default, as SAS also do. Can you write your model in mathematical terms, not all here understand STATA code. ... View more

hplogistic is in SAS/STAT, at least from  sas/stat 12.3. They are described in the documentation in its own section within the SAS/STAT. Correct that the Fisher scoring is used by default in proc logistic, but here the option "technique=newton" is Applied. ... View more

maybe, if you change to TECHNIQUE=NEWRAP in hplogistic, there is a change that it will use exactly same way from startingpoint approching the maximum. If not, then dont think more about it:-) ... View more

You should not worry about such small differences. The estimates calculated by most all procedures (including the proc logistic/hplogistic) is found by numerical maximation of a likelihood function. The difference you see  is caused by either different methods to maximize the function, or different tolerance of how close the solution should be to the maximum before it returns its results. If you are currious, you can try play a bit with the tolerence parameters. Try for instance to set the GCONV option at a smaller level.  Btw, when you set maxiter at 50, it just Means it allows for 50 iterations, but it will stop when it reach convergence criterium. Good luck. proc hplogistic data=&modelin. technique=NRRIDG GCONV=1E-6;   model &targetvar. (descending) = &modnum. / lackfit rsquare;   selection method=stepwise(slentry=&entry. slstay=&exit.) details=all;   output out=&modelout. p=pred copyvar=(&keyCol. &weightvar.); run; proc logistic data=&modelin. descending outest=lib.modparams;   weight &weightvar.;   model &targetvar.= &modnum.             /selection=stepwise slentry=&entry. slstay=&exit. RSQ lackfit outroc=Roc technique=newton MAXITER=50 GCONV=1E-6;   output out=&modelout. p=pred; run; ... View more