Calibration plot for large number of observations

[CraigSmith]

Hi all,

Hope you're well - wondering if you could help me with a query I had for calibration plots. 

Am validating  a logistic regression model in a large dataset (n=384,616). I can get the ROC curve to work but the calibration plot won't display.

Would anyone be able to modify my code (below) so the calibration plot will work?

Many thanks, 

Craig 

PROC LOGISTIC DATA=BIO3 PLOTS = ROC PLOTS(MAXPOINTS=NONE) 
 plots=calibration(CLM ShowObs);
MODEL HNC(EVENT='HNC') = PROB / ROCCI GOF OUTROC = ROC;
RUN; 

[PaigeMiller]

Please show us the ENTIRE log for this PROC LOGISTIC code — that's every single line in the log for PROC LOGISTIC, every single character in the log for PROC LOGISTIC.

[CraigSmith]

Hi Paige, 

Sorry was on a course w limited access to device.

Of course,

here is the  corresponding log 

---

1643
1644  PROC LOGISTIC DATA=BIO3 PLOTS = ROC PLOTS(MAXPOINTS=NONE)
1645   plots=calibration(CLM ShowObs);
1646  MODEL HNC(EVENT='HNC') = PROB / ROCCI GOF OUTROC = ROC;
1647  RUN;

NOTE: PROC LOGISTIC is modeling the probability that HNC='HNC'.
NOTE: Convergence criterion (GCONV=1E-8) satisfied.
WARNING: This graph has too many graphical elements. You may
         not be able to get any vector graphics output and in
         that case, you can set your output format to an image
         type.
WARNING: This graph has too many graphical elements. You may
         not be able to get any vector graphics output and in
         that case, you can set your output format to an image
         type.
NOTE: The number of observations of the LOESS plot (384616)
      exceeds the limit of 5000. Specify the LOESSMAXOBS option
      of the ODS GRAPHICS statement to override the limit.
NOTE: The number of observations of the LOESS plot (384616)
      exceeds the limit of 5000. Specify the LOESSMAXOBS option
      of the ODS GRAPHICS statement to override the limit.
NOTE: There were 384616 observations read from the data set
      WORK.BIO3.
NOTE: The data set WORK.ROC has 77689 observations and 8
      variables.
NOTE: PROCEDURE LOGISTIC used (Total process time):
      real time           4:43.01
      cpu time            3:42.45

---

[Ksharp]

I think there are must too many obs in your GRAPH to suppress calibration plot.

Or you could use PROC SGPLOT+ LOESS to get it.
Check @Rick_SAS 's blog:

https://blogs.sas.com/content/iml/2018/05/14/calibration-plots-in-sas.html

https://blogs.sas.com/content/iml/2018/05/16/decile-calibration-plots-sas.html

https://blogs.sas.com/content/iml/2019/02/20/easier-calibration-plot-sas.html

[CraigSmith]

Hi Ksharp,

I see  - thank you. Would it be possible for you to write an example with it please? Not quite sure how to do it. 

Many thanks,

C

[Ksharp]

/*
I think your obs is too big to use LOESS to get calibration plot.
Try "Decile calibration plots" ,
The following code is from :
https://blogs.sas.com/content/iml/2018/05/16/decile-calibration-plots-sas.html
*/

data have;
 set sashelp.heart sashelp.heart sashelp.heart sashelp.heart sashelp.heart
     sashelp.heart sashelp.heart sashelp.heart sashelp.heart sashelp.heart;
if status='Alive' then y=0;else y=1;
run;


/* Use PROC LOGISTIC and output the predicted probabilities.
   Intentionally MISSPECIFY the model as linear. */
proc logistic data=have noprint;
   model y(event='1') = weight height;
   output out=LogiOut predicted=PredProb; /* save predicted probabilities in data set */
run;
 
/* To construct the decile calibration plot, identify deciles of the predicted prob. */
proc rank data=LogiOut out=LogiDecile groups=10;
   var PredProb;
   ranks Decile;
run;
 
/* Then compute the mean predicted prob and the empirical proportions (and CI) for each decile */
proc means data=LogiDecile noprint;
   class Decile;
   types Decile;
   var y PredProb;
   output out=LogiDecileOut mean=yMean PredProbMean
          lclm=yLower uclm=yUpper;
run;
 
title "Calibration Plot for Misspecified Model";
title2 "True Model Is Quadratic; Fit Is Linear";
proc sgplot data=LogiDecileOut noautolegend aspect=1;
   lineparm x=0 y=0 slope=1 / lineattrs=(color=grey pattern=dash) clip;
   *loess x=PredProbMean y=yMean;  /* if you want a smoother based on deciles */
   series x=PredProbMean y=yMean;  /* if you to connect the deciles */
   scatter x=PredProbMean y=yMean / yerrorlower=yLower yerrorupper=yUpper;
   yaxis label="Observed Probability of Outcome";
   xaxis label="Predicted Probability of Outcome";
run;

[CraigSmith]

Hey Ksharp,

Thanks for this - this code is a bit too complex for me but I'll see if our statistician can get this running. 

Many thanks,

Craig