Regression Action Set

(View the complete code for this example.)

This section contains PROC CAS code.

Note: Input data must be accessible in your CAS session, either as a CAS table or as a transient-scope table. A CAS table has a two-level name: the first level is your CAS engine libref, and the second level is the table name. You refer to this table in the CAS procedure by specifying only the second level. For more information about two-level names, see Chapter 2, Shared Concepts. A transient-scope table is called directly from the action and exists in memory for the duration of the action. For more information about accessing data, see SAS Viya: System Programming Guide. For more information about PROC CAS and programming in CASL, see SAS Cloud Analytic Services: CASL Programmer’s Guide and SAS Cloud Analytic Services: CASL Reference.

The following DATA step creates the getStarted data table, which consists of 1,000 observations on a dichotomous response variable (y), a character variable (C), and 10 numeric variables (x1–x10).

proc cas;
dataStep.runCode /
    code="
   data getStarted;
      nTotalObs=1000;
      drop c2 eta pr i rew nTotalObs nObsPerThread nExtras;
      call streaminit(1);
      nObsPerThread = int(nTotalObs/_nthreads_);
      nExtras       = mod(nTotalObs,_nthreads_);
      if _threadid_ <= nExtras then nObsPerThread = nObsPerThread + 1;
      do i=1 to nObsPerThread;
         id = (_threadid_ - 1) * nObsPerThread + i;
         if _threadid_ > nExtras then id = id + nExtras;
         rew = rand('rewind', id);
         x1=round(rand('normal')*5+10,.1); x2=round(7*rand('uniform'))/7;
         x3=round(rand('normal')*1+2,.1);  x4=round(50*rand('uniform'));
         x5=round(100*rand('uniform'));    x6=round(rand('normal')*.8+1.5,.1);
         x7=10*round(10*rand('uniform'));  x8=round(10*rand('uniform'))/10;
         x9=round(rand('normal')*3+5,.1);  x10=round(rand('normal')*2+3,.1);
         c2=rand('uniform');
         if      (c2<.1) then C='A'; else if (c2<.2) then C='B';
         else if (c2<.3) then C='C'; else if (c2<.4) then C='D';
         else if (c2<.5) then C='E'; else if (c2<.6) then C='F';
         else if (c2<.7) then C='G'; else if (c2<.8) then C='H';
         else if (c2<.9) then C='I'; else                 C='J';
         eta=1-x2-x8;
         pr= exp(eta)/(1+exp(eta));
         y=(rand('uniform') > pr);
         output;
      end;
   run;
   ",
    single="no";
run;

By default, the CAS server executes the preceding DATA step in parallel in every thread on every node in your session. For example, if your DATA step creates 10 observations and executes in parallel on a session that has three worker nodes each of which uses 16 threads, then the resulting CAS table has 160 observations on each worker node, yielding a CAS table with a total of 480 observations. The preceding DATA step uses this technique to create a table with a specified number of observations in parallel on your CAS server by dividing the desired number of observations (nTotalObs) by the total number of threads (_nthreads_). If the number of observations cannot be exactly divided by the number of threads, the code adds one observation to each of a subset of the threads until all observations are accounted for.

The following statements use the logistic action in the regression action set to fit a logistic model to these data. The table parameter names the input data table to be analyzed. The class parameter specifies that the variable C is a classification variable; this parameter corresponds to the CLASS statement in the LOGSELECT procedure. The model parameter specifies the model: the depvars subparameter specifies the variable y as the response variable, and the effects subparameter lists the classification variable C and the 10 numeric explanatory variables, x1–x10. The itHist subparameter in the optimization parameter list produces a table that summarizes the steps of the optimization. The outputTables parameter outputs the fitted parameter estimates to the mycas.pe data table.

proc cas;
regression.logistic /
    class={"C"},
    model={depvar="y",
           effects={"C", "x1", "x2", "x3", "x4", "x5", "x6", "x7", "x8", "x9",
                    "x10"}},
    optimization={itHist="summary"},
    outputTables={names={parameterestimates="pe"}},
    table="getStarted";
run;

The output from this analysis is displayed by default, and is presented in Output 22.5.1 through Output 22.5.10.

Output 22.5.1 displays the "Model Information" table. The logistic action uses a Newton-Raphson algorithm to model a binary distribution with a logit link function for the variable y.

Output 22.5.2 displays the "Number of Observations" table. All 1,000 observations in the data table are used in the analysis.

The "Response Profile" table in Output 22.5.3 shows the breakdown of the response variable levels by frequency. By default for binary data, the logistic action models the probability of the event with the lower-ordered value in the "Response Profile" table, as indicated by the note that follows the table. (You can use the response-variable subparameters in the model parameter to choose which value of the response variable to model.) In this example, the values that are represented by y = 0 are modeled as the "successes" in the Bernoulli experiments.

The classification variable C has 10 unique formatted levels that are displayed in the "Class Level Information" table in Output 22.5.4.

The "Iteration History" table is shown in Output 22.5.5. The Newton-Raphson algorithm with ridging converged after three iterations, not counting the initial setup iteration.

Output 22.5.6 displays the final convergence status of the Newton-Raphson algorithm. The gconv relative convergence criterion is satisfied.

Output 22.5.7 displays the "Dimensions" table for this model. This table summarizes some important sizes of various model components. For example, it shows that the design matrix has 21 columns, which correspond to 1 column for the intercept, 10 columns for the effect associated with the classification variable C, and 1 column each for the continuous variables x1–x10. However, the rank of the crossproducts matrix is only 20. Because the classification variable C uses GLM parameterization and because the model contains an intercept, there is one singularity in the crossproducts matrix of the model. Consequently, only 20 parameters enter the optimization.

Output 22.5.8 shows the global test for the null hypothesis that all model effects jointly do not affect the probability of success of the binary response. The test is significant (p < 0.0001). One or more of the model effects thus significantly affect the probability of observing an event.

The "Fit Statistics" table is shown in Output 22.5.9. The –2 log likelihood at the converged estimates is 1,316.7. You can use this value to compare the model to nested model alternatives by means of a likelihood ratio test. To compare models that are not nested, you can use information criteria such as AIC (Akaike’s information criterion), AICC (Akaike’s bias-corrected information criterion), and SBC (Schwarz Bayesian information criterion). These criteria penalize the –2 log likelihood for the number of parameters. Because of the large number of parameters, the discrepancy between the –2 log likelihood and AIC (for example) is substantial in this case.

However, the "Parameter Estimates" table in Output 22.5.10 shows that many parameters have fairly large p-values, indicating that one or more of the model effects might not be necessary.

Finally, the procedure displays the table in Output 22.5.11, which shows the amount of time (in seconds) that the logistic action spent on different tasks in the analysis.

Output 22.5.10 shows that the full model has many nonsignificant effects. You can perform forward selection on these data by using the following statements. Effects that provide the best improvement for the default selection criterion, SBC, are added until no more effects can improve the selection criterion. The method subparameter in the selection parameter enables you to choose the selection method, and the details subparameter produces all tables that are related to model selection. The traceNames subparameter in the display parameter displays the full pathnames of the output tables (not shown here).

proc cas;
regression.logistic /
    class={"C"},
    display={traceNames="true"},
    model={depvar="y",
           effects={"C", "x1", "x2", "x3", "x4", "x5", "x6", "x7", "x8", "x9",
                    "x10"}},
    selection={details="all", method="forward"},
    table="getStarted";
run;

The model selection tables are shown in Output 22.5.12 through Output 22.5.14. Results from the selected model are shown in Output 22.5.15 and Output 22.5.16.

The "Selection Information" table in Output 22.5.12 summarizes the settings for the model selection. At each step, an effect is added to the model only if it produces the most significant improvement in the SBC. The forward selection stops three steps after the smallest SBC is obtained, or when all effects have been added to the model.

For each step of the selection process, the details=all subparameter displays the effects that are candidates for entering the model along with their selection criterion (SBC). Output 22.5.13 displays this table for the first step; the other steps are not shown here.

The details=all subparameter also displays the dimensions, global test, fit statistics, and parameter estimates at each step of the selection process; these details are not shown here.

When the selection procedure is complete, the "Selection Summary" table in Output 22.5.14 shows the effects that were added to the model and the value of their selection criterion (and the choose and stop criteria, if they are specified). Step 0 refers to the null model that contains only an intercept. In step 1, effect x8 made the most significant contribution to the model among the candidate effects, according to the SBC statistic. In step 2, x2 made the most significant contribution when an effect was added to a model that contains the intercept and x8. In the three subsequent steps, no effect could be added to the model that would reduce the SBC, so variable selection stopped because the stop horizon indicates that at most three steps beyond the minimum SBC value are to be taken (see Output 22.5.12) .

In Output 22.5.14, the "Selection Summary" table is followed by three small tables that summarize why the process stopped and which model is selected.

Output 22.5.15 displays information about the selected model. Notice that the –2 log-likelihood value in the "Fit Statistics" table is larger than the value for the full model in Output 22.5.9. This is expected because the selected model contains only a subset of the parameters. Because the selected model is more parsimonious than the full model, the discrepancy between the –2 log likelihood and the information criteria is less severe than previously noted.

The parameter estimates of the selected model are shown in Output 22.5.16. Notice that the effects are listed in the "Parameter Estimates" table in the order in which they were specified in the model parameter and not in the order in which they were added to the model.

You can construct the prediction equation for this model from the parameter estimates as follows. The estimated linear predictor for an observation is

The predicted probability that variable y takes the value 0 is

(View the complete code for this example.)

This section contains Lua code that creates the same getStarted table and performs the same analysis as shown in the CASL version of this example. For more information about coding in Lua, see Getting Started with SAS Viya for Lua and SAS Viya: System Programming Guide.

Because this example uses simulated data that are produced from DATA step code, you can submit the runCode action on your Lua client to create the simulated data in your CAS session. The following DATA step creates the getStarted data table, which consists of 1,000 observations on a dichotomous response variable (y), a character variable (C), and 10 numeric variables (x1–x10).

s:dataStep_runCode{single='no', code=[[
data getStarted;
   nTotalObs=1000;
   drop c2 eta pr i rew nTotalObs nObsPerThread nExtras;
   call streaminit(1);
   nObsPerThread = int(nTotalObs/_nthreads_);
   nExtras       = mod(nTotalObs,_nthreads_);
   if _threadid_ <= nExtras then nObsPerThread = nObsPerThread + 1;
   do i=1 to nObsPerThread;
      id = (_threadid_ - 1) * nObsPerThread + i;
      if _threadid_ > nExtras then id = id + nExtras;
      rew = rand('rewind', id);
      x1=round(rand('normal')*5+10,.1); x2=round(7*rand('uniform'))/7;
      x3=round(rand('normal')*1+2,.1);  x4=round(50*rand('uniform'));
      x5=round(100*rand('uniform'));    x6=round(rand('normal')*.8+1.5,.1);
      x7=10*round(10*rand('uniform'));  x8=round(10*rand('uniform'))/10;
      x9=round(rand('normal')*3+5,.1);  x10=round(rand('normal')*2+3,.1);
      c2=rand('uniform');
      if      (c2<.1) then C='A'; else if (c2<.2) then C='B';
      else if (c2<.3) then C='C'; else if (c2<.4) then C='D';
      else if (c2<.5) then C='E'; else if (c2<.6) then C='F';
      else if (c2<.7) then C='G'; else if (c2<.8) then C='H';
      else if (c2<.9) then C='I'; else                 C='J';
      eta=1-x2-x8;
      pr= exp(eta)/(1+exp(eta));
      y=(rand('uniform') > pr);
      output;
   end;
run;
]] }

By default, the CAS server executes the preceding DATA step in parallel in every thread on every node in your session. For example, if your DATA step creates 10 observations and executes in parallel on a session that has three worker nodes each of which uses 16 threads, then the resulting CAS table has 160 observations on each worker node, yielding a CAS table with a total of 480 observations. The preceding DATA step uses this technique to create a table with a specified number of observations in parallel on your CAS server by dividing the desired number of observations (nTotalObs) by the total number of threads (_nthreads_). If the number of observations cannot be exactly divided by the number of threads, the code adds one observation to each of a subset of the threads until all observations are accounted for.

The following statements use the logistic action in the regression action set to fit a logistic model to these data. The table parameter names the input data table to be analyzed. The class parameter specifies that the variable C is a classification variable; this parameter corresponds to the CLASS statement in the LOGSELECT procedure. The model parameter specifies the model: the depvars subparameter specifies the variable y as the response variable, and the effects subparameter lists the classification variable C and the 10 numeric explanatory variables, x1–x10. The itHist subparameter in the optimization parameter list produces a table that summarizes the steps of the optimization. The outputTables parameter outputs the fitted parameter estimates to the pe data table.

s:loadactionset{actionset="regression"}

s:regression_logistic{
   table='getStarted',
   class={'C'},
   model={depvar='y',
          effects={'C', 'x1', 'x2', 'x3', 'x4', 'x5',
                   'x6', 'x7', 'x8', 'x9', 'x10'}},
   optimization={itHist='summary'},
   outputTables={names={parameterestimates='pe'}} }

The output from this analysis is displayed by default, and is presented in Output 22.5.17 through Output 22.5.26.

Output 22.5.17 displays the "Model Information" table. The logistic action uses a Newton-Raphson algorithm to model a binary distribution with a logit link function for the variable y.

Output 22.5.18 displays the "Number of Observations" table. All 1,000 observations in the data table are used in the analysis.

The "Response Profile" table in Output 22.5.19 shows the breakdown of the response variable levels by frequency. By default for binary data, the logistic action models the probability of the event with the lower-ordered value in the "Response Profile" table, as indicated by the note that follows the table. (You can use the response-variable subparameters in the model parameter to choose which value of the response variable to model.) In this example, the values that are represented by y = 0 are modeled as the "successes" in the Bernoulli experiments.

The classification variable C has 10 unique formatted levels that are displayed in the "Class Level Information" table in Output 22.5.20.

The "Iteration History" table is shown in Output 22.5.21. The Newton-Raphson algorithm with ridging converged after three iterations, not counting the initial setup iteration.

Output 22.5.22 displays the final convergence status of the Newton-Raphson algorithm. The gconv relative convergence criterion is satisfied.

Output 22.5.23 displays the "Dimensions" table for this model. This table summarizes some important sizes of various model components. For example, it shows that the design matrix has 21 columns, which correspond to 1 column for the intercept, 10 columns for the effect associated with the classification variable C, and 1 column each for the continuous variables x1–x10. However, the rank of the crossproducts matrix is only 20. Because the classification variable C uses GLM parameterization and because the model contains an intercept, there is one singularity in the crossproducts matrix of the model. Consequently, only 20 parameters enter the optimization.

Output 22.5.24 shows the global test for the null hypothesis that all model effects jointly do not affect the probability of success of the binary response. The test is significant (p < 0.0001). One or more of the model effects thus significantly affect the probability of observing an event.

The "Fit Statistics" table is shown in Output 22.5.25. The –2 log likelihood at the converged estimates is 1,316.7. You can use this value to compare the model to nested model alternatives by means of a likelihood ratio test. To compare models that are not nested, you can use information criteria such as AIC (Akaike’s information criterion), AICC (Akaike’s bias-corrected information criterion), and SBC (Schwarz Bayesian information criterion). These criteria penalize the –2 log likelihood for the number of parameters. Because of the large number of parameters, the discrepancy between the –2 log likelihood and AIC (for example) is substantial in this case.

However, the "Parameter Estimates" table in Output 22.5.26 shows that many parameters have fairly large p-values, indicating that one or more of the model effects might not be necessary.

Finally, the procedure displays the table in Output 22.5.27, which shows the amount of time (in seconds) that the logistic action spent on different tasks in the analysis.

Output 22.5.26 shows that the full model has many nonsignificant effects. You can perform forward selection on these data by using the following statements. Effects that provide the best improvement for the default selection criterion, SBC, are added until no more effects can improve the selection criterion. The method subparameter in the selection parameter enables you to choose the selection method, and the details subparameter produces all tables that are related to model selection. The traceNames subparameter in the display parameter displays the full pathnames of the output tables (not shown here).

s:regression_logistic{
   table='getStarted',
   class={'C'},
   model={depvar='y',
          effects={'C', 'x1', 'x2', 'x3', 'x4', 'x5',
                   'x6', 'x7', 'x8', 'x9', 'x10'}},
   selection={method='forward', details='all'},
   display={traceNames=true} }

The model selection tables are shown in Output 22.5.28 through Output 22.5.30. Results from the selected model are shown in Output 22.5.31 and Output 22.5.32.

The "Selection Information" table in Output 22.5.28 summarizes the settings for the model selection. At each step, an effect is added to the model only if it produces the most significant improvement in the SBC. The forward selection stops three steps after the smallest SBC is obtained, or when all effects have been added to the model.

For each step of the selection process, the details=all subparameter displays the effects that are candidates for entering the model along with their selection criterion (SBC). Output 22.5.29 displays this table for the first step; the other steps are not shown here.

The details=all subparameter also displays the dimensions, global test, fit statistics, and parameter estimates at each step of the selection process; these details are not shown here.

When the selection procedure is complete, the "Selection Summary" table in Output 22.5.30 shows the effects that were added to the model and the value of their selection criterion (and the choose and stop criteria, if they are specified). Step 0 refers to the null model that contains only an intercept. In step 1, effect x8 made the most significant contribution to the model among the candidate effects, according to the SBC statistic. In step 2, x2 made the most significant contribution when an effect was added to a model that contains the intercept and x8. In the three subsequent steps, no effect could be added to the model that would reduce the SBC, so variable selection stopped because the stop horizon indicates that at most three steps beyond the minimum SBC value are to be taken (see Output 22.5.28) .

In Output 22.5.30, the "Selection Summary" table is followed by three small tables that summarize why the process stopped and which model is selected.

Output 22.5.31 displays information about the selected model. Notice that the –2 log-likelihood value in the "Fit Statistics" table is larger than the value for the full model in Output 22.5.25. This is expected because the selected model contains only a subset of the parameters. Because the selected model is more parsimonious than the full model, the discrepancy between the –2 log likelihood and the information criteria is less severe than previously noted.

The parameter estimates of the selected model are shown in Output 22.5.32. Notice that the effects are listed in the "Parameter Estimates" table in the order in which they were specified in the model parameter and not in the order in which they were added to the model.

You can construct the prediction equation for this model from the parameter estimates as follows. The estimated linear predictor for an observation is

The predicted probability that variable y takes the value 0 is

(View the complete code for this example.)

This section contains Python code that creates the same getStarted table and performs the same analysis as shown in the CASL version of this example. For more information about coding in Python, see Getting Started with SAS Viya for Python and SAS Viya: System Programming Guide.

Because this example uses simulated data that are produced from DATA step code, you can submit the runCode action on your Python client to create the simulated data in your CAS session. The following DATA step creates the getStarted data table, which consists of 1,000 observations on a dichotomous response variable (y), a character variable (C), and 10 numeric variables (x1–x10).

s.datastep.runcode(
    code='''
   data getStarted;
      nTotalObs=1000;
      drop c2 eta pr i rew nTotalObs nObsPerThread nExtras;
      call streaminit(1);
      nObsPerThread = int(nTotalObs/_nthreads_);
      nExtras       = mod(nTotalObs,_nthreads_);
      if _threadid_ <= nExtras then nObsPerThread = nObsPerThread + 1;
      do i=1 to nObsPerThread;
         id = (_threadid_ - 1) * nObsPerThread + i;
         if _threadid_ > nExtras then id = id + nExtras;
         rew = rand('rewind', id);
         x1=round(rand('normal')*5+10,.1); x2=round(7*rand('uniform'))/7;
         x3=round(rand('normal')*1+2,.1);  x4=round(50*rand('uniform'));
         x5=round(100*rand('uniform'));    x6=round(rand('normal')*.8+1.5,.1);
         x7=10*round(10*rand('uniform'));  x8=round(10*rand('uniform'))/10;
         x9=round(rand('normal')*3+5,.1);  x10=round(rand('normal')*2+3,.1);
         c2=rand('uniform');
         if      (c2<.1) then C='A'; else if (c2<.2) then C='B';
         else if (c2<.3) then C='C'; else if (c2<.4) then C='D';
         else if (c2<.5) then C='E'; else if (c2<.6) then C='F';
         else if (c2<.7) then C='G'; else if (c2<.8) then C='H';
         else if (c2<.9) then C='I'; else                 C='J';
         eta=1-x2-x8;
         pr= exp(eta)/(1+exp(eta));
         y=(rand('uniform') > pr);
         output;
      end;
   run;
   ''',
    single='no')

By default, the CAS server executes the preceding DATA step in parallel in every thread on every node in your session. For example, if your DATA step creates 10 observations and executes in parallel on a session that has three worker nodes each of which uses 16 threads, then the resulting CAS table has 160 observations on each worker node, yielding a CAS table with a total of 480 observations. The preceding DATA step uses this technique to create a table with a specified number of observations in parallel on your CAS server by dividing the desired number of observations (nTotalObs) by the total number of threads (_nthreads_). If the number of observations cannot be exactly divided by the number of threads, the code adds one observation to each of a subset of the threads until all observations are accounted for.

The following statements use the logistic action in the regression action set to fit a logistic model to these data. The table parameter names the input data table to be analyzed. The class parameter specifies that the variable C is a classification variable; this parameter corresponds to the CLASS statement in the LOGSELECT procedure. The model parameter specifies the model: the depvars subparameter specifies the variable y as the response variable, and the effects subparameter lists the classification variable C and the 10 numeric explanatory variables, x1–x10. The itHist subparameter in the optimization parameter list produces a table that summarizes the steps of the optimization. The outputTables parameter outputs the fitted parameter estimates to the pe data table.

s.loadactionset(actionset='regression')

s.regression.logistic(
    class_=['C'],
    model=dict(depvar='y',
               effects=['C', 'x1', 'x2', 'x3', 'x4', 'x5', 'x6', 'x7', 'x8',
                        'x9', 'x10']),
    optimization=dict(itHist='summary'),
    outputTables=dict(names=dict(parameterestimates='pe')),
    table='getStarted')

The output from this analysis is displayed by default, and is presented in Output 22.5.33 through Output 22.5.42.

Output 22.5.33 displays the "Model Information" table. The logistic action uses a Newton-Raphson algorithm to model a binary distribution with a logit link function for the variable y.

Output 22.5.34 displays the "Number of Observations" table. All 1,000 observations in the data table are used in the analysis.

The "Response Profile" table in Output 22.5.35 shows the breakdown of the response variable levels by frequency. By default for binary data, the logistic action models the probability of the event with the lower-ordered value in the "Response Profile" table, as indicated by the note that follows the table. (You can use the response-variable subparameters in the model parameter to choose which value of the response variable to model.) In this example, the values that are represented by y = 0 are modeled as the "successes" in the Bernoulli experiments.

The classification variable C has 10 unique formatted levels that are displayed in the "Class Level Information" table in Output 22.5.36.

The "Iteration History" table is shown in Output 22.5.37. The Newton-Raphson algorithm with ridging converged after three iterations, not counting the initial setup iteration.

Output 22.5.38 displays the final convergence status of the Newton-Raphson algorithm. The gconv relative convergence criterion is satisfied.

Output 22.5.39 displays the "Dimensions" table for this model. This table summarizes some important sizes of various model components. For example, it shows that the design matrix has 21 columns, which correspond to 1 column for the intercept, 10 columns for the effect associated with the classification variable C, and 1 column each for the continuous variables x1–x10. However, the rank of the crossproducts matrix is only 20. Because the classification variable C uses GLM parameterization and because the model contains an intercept, there is one singularity in the crossproducts matrix of the model. Consequently, only 20 parameters enter the optimization.

Output 22.5.40 shows the global test for the null hypothesis that all model effects jointly do not affect the probability of success of the binary response. The test is significant (p < 0.0001). One or more of the model effects thus significantly affect the probability of observing an event.

The "Fit Statistics" table is shown in Output 22.5.41. The –2 log likelihood at the converged estimates is 1,316.7. You can use this value to compare the model to nested model alternatives by means of a likelihood ratio test. To compare models that are not nested, you can use information criteria such as AIC (Akaike’s information criterion), AICC (Akaike’s bias-corrected information criterion), and SBC (Schwarz Bayesian information criterion). These criteria penalize the –2 log likelihood for the number of parameters. Because of the large number of parameters, the discrepancy between the –2 log likelihood and AIC (for example) is substantial in this case.

However, the "Parameter Estimates" table in Output 22.5.42 shows that many parameters have fairly large p-values, indicating that one or more of the model effects might not be necessary.

Finally, the procedure displays the table in Output 22.5.43, which shows the amount of time (in seconds) that the logistic action spent on different tasks in the analysis.

Output 22.5.42 shows that the full model has many nonsignificant effects. You can perform forward selection on these data by using the following statements. Effects that provide the best improvement for the default selection criterion, SBC, are added until no more effects can improve the selection criterion. The method subparameter in the selection parameter enables you to choose the selection method, and the details subparameter produces all tables that are related to model selection. The traceNames subparameter in the display parameter displays the full pathnames of the output tables (not shown here).

s.regression.logistic(
    class_=['C'],
    display=dict(traceNames='true'),
    model=dict(depvar='y',
               effects=['C', 'x1', 'x2', 'x3', 'x4', 'x5', 'x6', 'x7', 'x8',
                        'x9', 'x10']),
    selection=dict(details='all', method='forward'),
    table='getStarted')

The model selection tables are shown in Output 22.5.44 through Output 22.5.46. Results from the selected model are shown in Output 22.5.47 and Output 22.5.48.

The "Selection Information" table in Output 22.5.44 summarizes the settings for the model selection. At each step, an effect is added to the model only if it produces the most significant improvement in the SBC. The forward selection stops three steps after the smallest SBC is obtained, or when all effects have been added to the model.

For each step of the selection process, the details=all subparameter displays the effects that are candidates for entering the model along with their selection criterion (SBC). Output 22.5.45 displays this table for the first step; the other steps are not shown here.

The details=all subparameter also displays the dimensions, global test, fit statistics, and parameter estimates at each step of the selection process; these details are not shown here.

When the selection procedure is complete, the "Selection Summary" table in Output 22.5.46 shows the effects that were added to the model and the value of their selection criterion (and the choose and stop criteria, if they are specified). Step 0 refers to the null model that contains only an intercept. In step 1, effect x8 made the most significant contribution to the model among the candidate effects, according to the SBC statistic. In step 2, x2 made the most significant contribution when an effect was added to a model that contains the intercept and x8. In the three subsequent steps, no effect could be added to the model that would reduce the SBC, so variable selection stopped because the stop horizon indicates that at most three steps beyond the minimum SBC value are to be taken (see Output 22.5.44) .

In Output 22.5.46, the "Selection Summary" table is followed by three small tables that summarize why the process stopped and which model is selected.

Output 22.5.47 displays information about the selected model. Notice that the –2 log-likelihood value in the "Fit Statistics" table is larger than the value for the full model in Output 22.5.41. This is expected because the selected model contains only a subset of the parameters. Because the selected model is more parsimonious than the full model, the discrepancy between the –2 log likelihood and the information criteria is less severe than previously noted.

The parameter estimates of the selected model are shown in Output 22.5.48. Notice that the effects are listed in the "Parameter Estimates" table in the order in which they were specified in the model parameter and not in the order in which they were added to the model.

You can construct the prediction equation for this model from the parameter estimates as follows. The estimated linear predictor for an observation is

The predicted probability that variable y takes the value 0 is

(View the complete code for this example.)

This section contains R code that creates the same getStarted table and performs the same analysis as shown in the CASL version of this example. For more information about coding in R, see Getting Started with SAS Viya for R and SAS Viya: System Programming Guide.

Because this example uses simulated data that are produced from DATA step code, you can submit the runCode action on your R client to create the simulated data in your CAS session. The following DATA step creates the getStarted data table, which consists of 1,000 observations on a dichotomous response variable (y), a character variable (C), and 10 numeric variables (x1–x10).

cas.dataStep.runCode(s,
    code='
   data getStarted;
      nTotalObs=1000;
      drop c2 eta pr i rew nTotalObs nObsPerThread nExtras;
      call streaminit(1);
      nObsPerThread = int(nTotalObs/_nthreads_);
      nExtras       = mod(nTotalObs,_nthreads_);
      if _threadid_ <= nExtras then nObsPerThread = nObsPerThread + 1;
      do i=1 to nObsPerThread;
         id = (_threadid_ - 1) * nObsPerThread + i;
         if _threadid_ > nExtras then id = id + nExtras;
         rew = rand(\'rewind\', id);
         x1=round(rand(\'normal\')*5+10,.1); x2=round(7*rand(\'uniform\'))/7;
         x3=round(rand(\'normal\')*1+2,.1);  x4=round(50*rand(\'uniform\'));
         x5=round(100*rand(\'uniform\'));
      x6=round(rand(\'normal\')*.8+1.5,.1);
         x7=10*round(10*rand(\'uniform\'));  x8=round(10*rand(\'uniform\'))/10;
         x9=round(rand(\'normal\')*3+5,.1);  x10=round(rand(\'normal\')*2+3,.1);
         c2=rand(\'uniform\');
         if      (c2<.1) then C=\'A\'; else if (c2<.2) then C=\'B\';
         else if (c2<.3) then C=\'C\'; else if (c2<.4) then C=\'D\';
         else if (c2<.5) then C=\'E\'; else if (c2<.6) then C=\'F\';
         else if (c2<.7) then C=\'G\'; else if (c2<.8) then C=\'H\';
         else if (c2<.9) then C=\'I\'; else                 C=\'J\';
         eta=1-x2-x8;
         pr= exp(eta)/(1+exp(eta));
         y=(rand(\'uniform\') > pr);
         output;
      end;
   run;
   ',
    single='no')

By default, the CAS server executes the preceding DATA step in parallel in every thread on every node in your session. For example, if your DATA step creates 10 observations and executes in parallel on a session that has three worker nodes each of which uses 16 threads, then the resulting CAS table has 160 observations on each worker node, yielding a CAS table with a total of 480 observations. The preceding DATA step uses this technique to create a table with a specified number of observations in parallel on your CAS server by dividing the desired number of observations (nTotalObs) by the total number of threads (_nthreads_). If the number of observations cannot be exactly divided by the number of threads, the code adds one observation to each of a subset of the threads until all observations are accounted for.

The following statements use the logistic action in the regression action set to fit a logistic model to these data. The table parameter names the input data table to be analyzed. The class parameter specifies that the variable C is a classification variable; this parameter corresponds to the CLASS statement in the LOGSELECT procedure. The model parameter specifies the model: the depvars subparameter specifies the variable y as the response variable, and the effects subparameter lists the classification variable C and the 10 numeric explanatory variables, x1–x10. The itHist subparameter in the optimization parameter list produces a table that summarizes the steps of the optimization. The outputTables parameter outputs the fitted parameter estimates to the pe data table.

m <- cas.builtins.loadActionSet(s, actionset='regression')

cas.regression.logistic(s,
    class=c('C'),
    model=list(depvar='y',
               effects=c('C', 'x1', 'x2', 'x3', 'x4', 'x5', 'x6', 'x7', 'x8',
                         'x9', 'x10')),
    optimization=list(itHist='summary'),
    outputTables=list(names=list(parameterestimates='pe')),
    table='getStarted')

The output from this analysis is displayed by default, and is presented in Output 22.5.49 through Output 22.5.58.

Output 22.5.49 displays the "Model Information" table. The logistic action uses a Newton-Raphson algorithm to model a binary distribution with a logit link function for the variable y.

Output 22.5.50 displays the "Number of Observations" table. All 1,000 observations in the data table are used in the analysis.

The "Response Profile" table in Output 22.5.51 shows the breakdown of the response variable levels by frequency. By default for binary data, the logistic action models the probability of the event with the lower-ordered value in the "Response Profile" table, as indicated by the note that follows the table. (You can use the response-variable subparameters in the model parameter to choose which value of the response variable to model.) In this example, the values that are represented by y = 0 are modeled as the "successes" in the Bernoulli experiments.

The classification variable C has 10 unique formatted levels that are displayed in the "Class Level Information" table in Output 22.5.52.

The "Iteration History" table is shown in Output 22.5.53. The Newton-Raphson algorithm with ridging converged after three iterations, not counting the initial setup iteration.

Output 22.5.54 displays the final convergence status of the Newton-Raphson algorithm. The gconv relative convergence criterion is satisfied.

Output 22.5.55 displays the "Dimensions" table for this model. This table summarizes some important sizes of various model components. For example, it shows that the design matrix has 21 columns, which correspond to 1 column for the intercept, 10 columns for the effect associated with the classification variable C, and 1 column each for the continuous variables x1–x10. However, the rank of the crossproducts matrix is only 20. Because the classification variable C uses GLM parameterization and because the model contains an intercept, there is one singularity in the crossproducts matrix of the model. Consequently, only 20 parameters enter the optimization.

Output 22.5.56 shows the global test for the null hypothesis that all model effects jointly do not affect the probability of success of the binary response. The test is significant (p < 0.0001). One or more of the model effects thus significantly affect the probability of observing an event.

The "Fit Statistics" table is shown in Output 22.5.57. The –2 log likelihood at the converged estimates is 1,316.7. You can use this value to compare the model to nested model alternatives by means of a likelihood ratio test. To compare models that are not nested, you can use information criteria such as AIC (Akaike’s information criterion), AICC (Akaike’s bias-corrected information criterion), and SBC (Schwarz Bayesian information criterion). These criteria penalize the –2 log likelihood for the number of parameters. Because of the large number of parameters, the discrepancy between the –2 log likelihood and AIC (for example) is substantial in this case.

However, the "Parameter Estimates" table in Output 22.5.58 shows that many parameters have fairly large p-values, indicating that one or more of the model effects might not be necessary.

Finally, the procedure displays the table in Output 22.5.59, which shows the amount of time (in seconds) that the logistic action spent on different tasks in the analysis.

Output 22.5.58 shows that the full model has many nonsignificant effects. You can perform forward selection on these data by using the following statements. Effects that provide the best improvement for the default selection criterion, SBC, are added until no more effects can improve the selection criterion. The method subparameter in the selection parameter enables you to choose the selection method, and the details subparameter produces all tables that are related to model selection. The traceNames subparameter in the display parameter displays the full pathnames of the output tables (not shown here).

cas.regression.logistic(s,
    class=c('C'),
    display=list(traceNames='true'),
    model=list(depvar='y',
               effects=c('C', 'x1', 'x2', 'x3', 'x4', 'x5', 'x6', 'x7', 'x8',
                         'x9', 'x10')),
    selection=list(details='all', method='forward'),
    table='getStarted')

The model selection tables are shown in Output 22.5.60 through Output 22.5.62. Results from the selected model are shown in Output 22.5.63 and Output 22.5.64.

The "Selection Information" table in Output 22.5.60 summarizes the settings for the model selection. At each step, an effect is added to the model only if it produces the most significant improvement in the SBC. The forward selection stops three steps after the smallest SBC is obtained, or when all effects have been added to the model.

For each step of the selection process, the details=all subparameter displays the effects that are candidates for entering the model along with their selection criterion (SBC). Output 22.5.61 displays this table for the first step; the other steps are not shown here.

The details=all subparameter also displays the dimensions, global test, fit statistics, and parameter estimates at each step of the selection process; these details are not shown here.

When the selection procedure is complete, the "Selection Summary" table in Output 22.5.62 shows the effects that were added to the model and the value of their selection criterion (and the choose and stop criteria, if they are specified). Step 0 refers to the null model that contains only an intercept. In step 1, effect x8 made the most significant contribution to the model among the candidate effects, according to the SBC statistic. In step 2, x2 made the most significant contribution when an effect was added to a model that contains the intercept and x8. In the three subsequent steps, no effect could be added to the model that would reduce the SBC, so variable selection stopped because the stop horizon indicates that at most three steps beyond the minimum SBC value are to be taken (see Output 22.5.60) .

In Output 22.5.62, the "Selection Summary" table is followed by three small tables that summarize why the process stopped and which model is selected.

Output 22.5.63 displays information about the selected model. Notice that the –2 log-likelihood value in the "Fit Statistics" table is larger than the value for the full model in Output 22.5.57. This is expected because the selected model contains only a subset of the parameters. Because the selected model is more parsimonious than the full model, the discrepancy between the –2 log likelihood and the information criteria is less severe than previously noted.

The parameter estimates of the selected model are shown in Output 22.5.64. Notice that the effects are listed in the "Parameter Estimates" table in the order in which they were specified in the model parameter and not in the order in which they were added to the model.

You can construct the prediction equation for this model from the parameter estimates as follows. The estimated linear predictor for an observation is

The predicted probability that variable y takes the value 0 is

(View the complete code for this example.)

This section contains PROC CAS code.

Note: Input data must be accessible in your CAS session, either as a CAS table or as a transient-scope table. A CAS table has a two-level name: the first level is your CAS engine libref, and the second level is the table name. You refer to this table in the CAS procedure by specifying only the second level. For more information about two-level names, see Chapter 2, Shared Concepts. A transient-scope table is called directly from the action and exists in memory for the duration of the action. For more information about accessing data, see SAS Viya: System Programming Guide. For more information about PROC CAS and programming in CASL, see SAS Cloud Analytic Services: CASL Programmer’s Guide and SAS Cloud Analytic Services: CASL Reference.

The following DATA step creates the data table mycas.Ingots—which consists of the number, r, of ingots not ready for rolling, out of n tested, for a number of combinations of heating time (Heat) and soaking time (Soak)—in your CAS session. This DATA step assumes that your CAS engine libref is named mycas, but you can substitute any appropriately defined CAS engine libref.

data mycas.Ingots;
   input Heat Soak r n @@;
   a = n - r;
   Obsnum = _n_;
   datalines;
7 1.0 0 10  14 1.0 0 31  27 1.0 1 56  51 1.0 3 13
7 1.7 0 17  14 1.7 0 43  27 1.7 4 44  51 1.7 0  1
7 2.2 0  7  14 2.2 2 33  27 2.2 0 21  51 2.2 0  1
7 2.8 0 12  14 2.8 0 31  27 2.8 1 22  51 4.0 0  1
7 4.0 0  9  14 4.0 0 19  27 4.0 1 16
;

The following statements use events/trials syntax to model the binomial response. The events variable in this situation is r, and the trials variable is n. The dependency of the probability of not being ready for rolling is modeled as a function of heating time (Heat), soaking time (Soak), and their interaction. The association parameter displays ordinal measures of association between the observed responses and the predicted probabilities. The ctable parameter saves statistics to the mycas.Roc data table that are used to evaluate the predictive power of the model. The output parameter stores the linear predictors and the predicted probabilities in the mycas.Out data table along with the other variables from the input data table.

proc cas;
regression.logistic /
    association="true",
    ctable={casOut={name="Roc", replace="true"},
            fpf="FPF",
            nocounts="true",
            tpf="TPF"},
    model={depVars={{name="r"}},
           effects={{vars={"Heat", "Soak"}},
                    {interaction="CROSS", vars={"Heat", "Soak"}}},
           trial="n"},
    output={casOut={name="Out", replace="true"},
            copyVars={"Heat", "Soak"},
            pred="Pred",
            xBeta="_XBETA_"},
    table="Ingots";
run;

The "Model Information" table shows that the data are modeled as binomially distributed with a logit link function (Output 22.6.1). This is the default link function in the logistic action for binary and binomial data. The action estimates the parameters of the model by a Newton-Raphson algorithm.

The second table in Output 22.6.1 shows that all 19 observations in the input data table were used in the analysis and the "Response Profile" table shows that the total number of events and nonevents equals 12 and 375, respectively.

Output 22.6.2 displays the convergence status table for this run. The logistic action satisfies the gconv convergence criterion.

Output 22.6.3 displays the "Dimensions" table for the model. The design matrix of the model (the matrix) has four columns, which correspond to the intercept, the Heat effect, the Soak effect, and the interaction of the Heat and Soak effects. The model is nonsingular because the rank of the crossproducts matrix equals the number of columns in . All parameters are estimable, and all participate in the optimization.

Output 22.6.4 displays the "Fit Statistics" table for this run. Evaluated at the converged estimates, –2 times the value of the log-likelihood function equals 27.95689. Further fit statistics are also displayed, all of them in "smaller is better" form. The AIC, AICC, and SBC criteria are used to compare non-nested models and to penalize the model fit for the number of observations and parameters. You can use the –2 log-likelihood value to compare nested models by way of a likelihood ratio test.

Output 22.6.5 shows the test of the global hypothesis that the effects jointly do not affect the probability of ingot readiness. You can obtain the chi-square test statistic by comparing the –2 log-likelihood value of the model with covariates to the value in the intercept-only model. The test is significant, with a p-value of 0.0082. One or more of the effects in the model have a significant impact on the probability of ingot readiness.

The "Parameter Estimates" table in Output 22.6.6 displays the estimates and standard errors of the model effects.

Output 22.6.7 displays the "Association Statistics" table, which is produced when you specify the association parameter. The table contains four measures of association for assessing the predictive ability of a model. For more information, see the section Association Statistics.

You can construct the prediction equation of the model from the "Parameter Estimates" table. For example, an observation with Heat equal to 14 and Soak equal to 1.7 has the linear predictor

The probability that an ingot with these characteristics is not ready for rolling is

The logistic action computes these linear predictors and probabilities and stores them in the mycas.Out data table.

The following fetch action from the tables action set enables you to display Output 22.6.8, which shows the probability that an ingot with Heat equal to 14 and Soak equal to 1.7 is not ready for rolling:

proc cas;
table.fetch /
    format="true",
    index="false",
    table={name="Out", where="14 = Heat and 1.7 = Soak"};
run;

The ctable parameter computes statistics for binary response models (on the basis of classifying observations according to whether their predicted probabilities exceed certain values) and stores the results in the mycas.Roc data table. For more information, see the section Classification Table and ROC Curves.

You can use the mycas.Roc data table to display the ROC curve by using the SGPLOT procedure as follows:

ods graphics on;
proc sgplot data=mycas.Roc aspect=1 noautolegend;
   title 'ROC Curve';
   xaxis values=(0 to 1 by 0.25) grid offsetmin=.05 offsetmax=.05;
   yaxis values=(0 to 1 by 0.25) grid offsetmin=.05 offsetmax=.05;
   lineparm x=0 y=0 slope=1 / lineattrs=(color=ligr);
   series x=FPF y=TPF;
   inset 'Area under the curve=0.7706' / position=bottomright;
run;

If you have counts of events and nonevents, then instead of using events/trials syntax, you can use the multinomial-trial syntax to perform the same analysis. In the mycas.Ingots data table, the a variable contains the number of ingots ready for rolling (the nonevents). In the following program, you specify the event variable, followed by the nonevent variable:

proc cas;
regression.logistic /
    model={depVars={{name="r"}, {name="a"}},
           effects={{vars={"Heat", "Soak"}},
                    {interaction="CROSS", vars={"Heat", "Soak"}}}},
    table="Ingots";
run;

The only difference from the preceding results is that in this case the "Model Information" table displays the event and nonevent variables; this table is not displayed here.

Binomial data are a form of grouped binary data in which "successes" in the underlying Bernoulli trials are totaled. You can thus unwind data for which you use the events/trials syntax and fit these data by using techniques for binary data.

The following DATA step expands the mycas.Ingots data table, which contains 12 events in 387 trials, into a data table that has 387 observations:

proc cas;
dataStep.runCode /
    code="
   data Ingots_binary;
      set Ingots;
      do i=1 to n;
        if i <= r then y=1; else y = 0;
        output;
      end;
   run;
   ",
    single="no";
run;

The following statements fit the model that contains the Heat effect, the Soak effect, and their interaction to the data table. The event=’1’ response-variable subparameter in the model parameter ensures that the logistic action models the probability that the variable y takes the value 1.

proc cas;
regression.logistic /
    model={depVars={{name="y", options={event="1"}}},
           effects={{vars={"Heat", "Soak"}},
                    {interaction="CROSS", vars={"Heat", "Soak"}}}},
    table={name="Ingots_binary"};
run;

Output 22.6.10 displays the "Model Information," "Number of Observations," and "Response Profile" tables. The data are now modeled as binary (Bernoulli distributed) with a logit link function. The "Response Profile" table shows that the binary response breaks down into 375 observations where y=0 and 12 observations where y=1.

Output 22.6.11 displays the result of the test of the global null hypothesis and the parameter estimates. These results match those in Output 22.6.5 and Output 22.6.6.

(View the complete code for this example.)

This section contains Lua code that performs the same analysis as shown in the CASL version of this example.

Note: In order to run this code, the data that are described in the CASL version need to be accessible to the CAS server. One way to do this is to convert the Ingots data to the comma-separated-value (CSV) file Ingots.csv and then use the following code to load the CSV file into CAS:

s:loadtable{casLib="casuser", path="Ingots.csv"}

For more information about coding in Lua, see Getting Started with SAS Viya for Lua and SAS Viya: System Programming Guide.

The following statements use events/trials syntax to model the binomial response. The events variable in this situation is r, and the trials variable is n. The dependency of the probability of not being ready for rolling is modeled as a function of heating time (Heat), soaking time (Soak), and their interaction. The association parameter displays ordinal measures of association between the observed responses and the predicted probabilities. The ctable parameter saves statistics to the Roc data table that are used to evaluate the predictive power of the model. The output parameter stores the linear predictors and the predicted probabilities in the Out data table along with the other variables from the input data table.

s:loadactionset{actionset="regression"}

s:regression.logistic {
   table='Ingots',
   model={depVars={{name='r'}},
          effects={{vars={'Heat', 'Soak'}},
                   {vars={'Heat', 'Soak'}, interaction='CROSS'}},
                    trial='n'},
   output={casOut={name='Out', replace=true},
           copyVars={'Heat', 'Soak'},
           pred='Pred', xBeta='_XBETA_'},
   association=true,
   ctable={casOut={name='Roc', replace=true},
           nocounts=true, tpf='TPF', fpf='FPF'} }

The "Model Information" table shows that the data are modeled as binomially distributed with a logit link function (Output 22.6.12). This is the default link function in the logistic action for binary and binomial data. The action estimates the parameters of the model by a Newton-Raphson algorithm.

The second table in Output 22.6.12 shows that all 19 observations in the input data table were used in the analysis and the "Response Profile" table shows that the total number of events and nonevents equals 12 and 375, respectively.

Output 22.6.13 displays the convergence status table for this run. The logistic action satisfies the gconv convergence criterion.

Output 22.6.14 displays the "Dimensions" table for the model. The design matrix of the model (the matrix) has four columns, which correspond to the intercept, the Heat effect, the Soak effect, and the interaction of the Heat and Soak effects. The model is nonsingular because the rank of the crossproducts matrix equals the number of columns in . All parameters are estimable, and all participate in the optimization.

Output 22.6.15 displays the "Fit Statistics" table for this run. Evaluated at the converged estimates, –2 times the value of the log-likelihood function equals 27.95689. Further fit statistics are also displayed, all of them in "smaller is better" form. The AIC, AICC, and SBC criteria are used to compare non-nested models and to penalize the model fit for the number of observations and parameters. You can use the –2 log-likelihood value to compare nested models by way of a likelihood ratio test.

Output 22.6.16 shows the test of the global hypothesis that the effects jointly do not affect the probability of ingot readiness. You can obtain the chi-square test statistic by comparing the –2 log-likelihood value of the model with covariates to the value in the intercept-only model. The test is significant, with a p-value of 0.0082. One or more of the effects in the model have a significant impact on the probability of ingot readiness.

The "Parameter Estimates" table in Output 22.6.17 displays the estimates and standard errors of the model effects.

Output 22.6.18 displays the "Association Statistics" table, which is produced when you specify the association parameter. The table contains four measures of association for assessing the predictive ability of a model. For more information, see the section Association Statistics.

You can construct the prediction equation of the model from the "Parameter Estimates" table. For example, an observation with Heat equal to 14 and Soak equal to 1.7 has the linear predictor

The probability that an ingot with these characteristics is not ready for rolling is

The logistic action computes these linear predictors and probabilities and stores them in the Out data table.

The following fetch action from the tables action set enables you to display Output 22.6.19, which shows the probability that an ingot with Heat equal to 14 and Soak equal to 1.7 is not ready for rolling:

s:table_fetch{table={name='Out', where='14 = Heat and 1.7 = Soak'},
              index=false, format=true}

The ctable parameter computes statistics for binary response models (on the basis of classifying observations according to whether their predicted probabilities exceed certain values) and stores the results in the Roc data table. For more information, see the section Classification Table and ROC Curves.

The Roc data table contains information you can use to produce an ROC curve like the following:

If you have counts of events and nonevents, then instead of using events/trials syntax, you can use the multinomial-trial syntax to perform the same analysis. In the mycas.Ingots data table, the a variable contains the number of ingots ready for rolling (the nonevents). In the following program, you specify the event variable, followed by the nonevent variable:

s:regression.logistic {
   table='Ingots',
   model={depVars={{name='r'},{name='a'}},
          effects={{vars={'Heat', 'Soak'}},
                   {vars={'Heat', 'Soak'}, interaction='CROSS'}}} }

The only difference from the preceding results is that in this case the "Model Information" table displays the event and nonevent variables; this table is not displayed here.

Binomial data are a form of grouped binary data in which "successes" in the underlying Bernoulli trials are totaled. You can thus unwind data for which you use the events/trials syntax and fit these data by using techniques for binary data.

The following DATA step expands the Ingots data table, which contains 12 events in 387 trials, into a data table that has 387 observations:

s:dataStep_runCode{single='no', code=[[
data Ingots_binary;
   set Ingots;
   do i=1 to n;
     if i <= r then y=1; else y = 0;
     output;
   end;
run;
]] }

The following statements fit the model that contains the Heat effect, the Soak effect, and their interaction to the data table. The event=’1’ response-variable subparameter in the model parameter ensures that the logistic action models the probability that the variable y takes the value 1.

s:regression_logistic {
   table={name='Ingots_binary'},
   model={depVars={{name='y', options={event='1'}}},
          effects={{vars={'Heat', 'Soak'}},
                   {vars={'Heat', 'Soak'}, interaction='CROSS'}}}}

Output 22.6.21 displays the "Model Information," "Number of Observations," and "Response Profile" tables. The data are now modeled as binary (Bernoulli distributed) with a logit link function. The "Response Profile" table shows that the binary response breaks down into 375 observations where y=0 and 12 observations where y=1.

Output 22.6.22 displays the result of the test of the global null hypothesis and the parameter estimates. These results match those in Output 22.6.16 and Output 22.6.17.

(View the complete code for this example.)

This section contains Python code that performs the same analysis as shown in the CASL version of this example.

Note: In order to run this code, the data that are described in the CASL version need to be accessible to the CAS server. One way to do this is to convert the Ingots data to the comma-separated-value (CSV) file Ingots.csv and then use the following code to load the CSV file into CAS:

s.upload_file('Ingots.csv')

For more information about coding in Python, see Getting Started with SAS Viya for Python and SAS Viya: System Programming Guide.

The following statements use events/trials syntax to model the binomial response. The events variable in this situation is r, and the trials variable is n. The dependency of the probability of not being ready for rolling is modeled as a function of heating time (Heat), soaking time (Soak), and their interaction. The association parameter displays ordinal measures of association between the observed responses and the predicted probabilities. The ctable parameter saves statistics to the Roc data table that are used to evaluate the predictive power of the model. The output parameter stores the linear predictors and the predicted probabilities in the Out data table along with the other variables from the input data table.

s.loadactionset(actionset='regression')

s.regression.logistic(
    association='true',
    ctable=dict(casOut=dict(name='Roc', replace='true'),
                fpf='FPF',
                nocounts='true',
                tpf='TPF'),
    model=dict(depVars=[dict(name='r')],
               effects=[dict(vars=['Heat', 'Soak']),
                        dict(interaction='CROSS', vars=['Heat', 'Soak'])],
               trial='n'),
    output=dict(casOut=dict(name='Out', replace='true'),
                copyVars=['Heat', 'Soak'],
                pred='Pred',
                xBeta='_XBETA_'),
    table='Ingots')

The "Model Information" table shows that the data are modeled as binomially distributed with a logit link function (Output 22.6.23). This is the default link function in the logistic action for binary and binomial data. The action estimates the parameters of the model by a Newton-Raphson algorithm.

The second table in Output 22.6.23 shows that all 19 observations in the input data table were used in the analysis and the "Response Profile" table shows that the total number of events and nonevents equals 12 and 375, respectively.

Output 22.6.24 displays the convergence status table for this run. The logistic action satisfies the gconv convergence criterion.

Output 22.6.25 displays the "Dimensions" table for the model. The design matrix of the model (the matrix) has four columns, which correspond to the intercept, the Heat effect, the Soak effect, and the interaction of the Heat and Soak effects. The model is nonsingular because the rank of the crossproducts matrix equals the number of columns in . All parameters are estimable, and all participate in the optimization.

Output 22.6.26 displays the "Fit Statistics" table for this run. Evaluated at the converged estimates, –2 times the value of the log-likelihood function equals 27.95689. Further fit statistics are also displayed, all of them in "smaller is better" form. The AIC, AICC, and SBC criteria are used to compare non-nested models and to penalize the model fit for the number of observations and parameters. You can use the –2 log-likelihood value to compare nested models by way of a likelihood ratio test.

Output 22.6.27 shows the test of the global hypothesis that the effects jointly do not affect the probability of ingot readiness. You can obtain the chi-square test statistic by comparing the –2 log-likelihood value of the model with covariates to the value in the intercept-only model. The test is significant, with a p-value of 0.0082. One or more of the effects in the model have a significant impact on the probability of ingot readiness.

The "Parameter Estimates" table in Output 22.6.28 displays the estimates and standard errors of the model effects.

Output 22.6.29 displays the "Association Statistics" table, which is produced when you specify the association parameter. The table contains four measures of association for assessing the predictive ability of a model. For more information, see the section Association Statistics.

You can construct the prediction equation of the model from the "Parameter Estimates" table. For example, an observation with Heat equal to 14 and Soak equal to 1.7 has the linear predictor

The probability that an ingot with these characteristics is not ready for rolling is

The logistic action computes these linear predictors and probabilities and stores them in the Out data table.

The following fetch action from the tables action set enables you to display Output 22.6.30, which shows the probability that an ingot with Heat equal to 14 and Soak equal to 1.7 is not ready for rolling:

s.table.fetch(
    format='true',
    index='false',
    table=dict(name='Out', where='14 = Heat and 1.7 = Soak'))

The ctable parameter computes statistics for binary response models (on the basis of classifying observations according to whether their predicted probabilities exceed certain values) and stores the results in the Roc data table. For more information, see the section Classification Table and ROC Curves.

The Roc data table contains information you can use to produce an ROC curve like the following:

If you have counts of events and nonevents, then instead of using events/trials syntax, you can use the multinomial-trial syntax to perform the same analysis. In the mycas.Ingots data table, the a variable contains the number of ingots ready for rolling (the nonevents). In the following program, you specify the event variable, followed by the nonevent variable:

s.regression.logistic(
    model=dict(depVars=[dict(name='r'), dict(name='a')],
               effects=[dict(vars=['Heat', 'Soak']),
                        dict(interaction='CROSS', vars=['Heat', 'Soak'])]),
    table='Ingots')

The only difference from the preceding results is that in this case the "Model Information" table displays the event and nonevent variables; this table is not displayed here.

Binomial data are a form of grouped binary data in which "successes" in the underlying Bernoulli trials are totaled. You can thus unwind data for which you use the events/trials syntax and fit these data by using techniques for binary data.

The following DATA step expands the Ingots data table, which contains 12 events in 387 trials, into a data table that has 387 observations:

s.datastep.runcode(
    code='''
   data Ingots_binary;
      set Ingots;
      do i=1 to n;
        if i <= r then y=1; else y = 0;
        output;
      end;
   run;
   ''',
    single='no')

The following statements fit the model that contains the Heat effect, the Soak effect, and their interaction to the data table. The event=’1’ response-variable subparameter in the model parameter ensures that the logistic action models the probability that the variable y takes the value 1.

s.regression.logistic(
    model=dict(depVars=[dict(name='y', options=dict(event='1'))],
               effects=[dict(vars=['Heat', 'Soak']),
                        dict(interaction='CROSS', vars=['Heat', 'Soak'])]),
    table=dict(name='Ingots_binary'))

Output 22.6.32 displays the "Model Information," "Number of Observations," and "Response Profile" tables. The data are now modeled as binary (Bernoulli distributed) with a logit link function. The "Response Profile" table shows that the binary response breaks down into 375 observations where y=0 and 12 observations where y=1.

Output 22.6.33 displays the result of the test of the global null hypothesis and the parameter estimates. These results match those in Output 22.6.27 and Output 22.6.28.

(View the complete code for this example.)

This section contains R code that performs the same analysis as shown in the CASL version of this example.

Note: In order to run this code, the data that are described in the CASL version need to be accessible to the CAS server. One way to do this is to convert the Ingots data to the comma-separated-value (CSV) file Ingots.csv and then use the following code to load the CSV file into CAS:

m <- cas.read.csv(s, "Ingots.csv", casOut=list(name="Ingots"))

For more information about coding in R, see Getting Started with SAS Viya for R and SAS Viya: System Programming Guide.

The following statements use events/trials syntax to model the binomial response. The events variable in this situation is r, and the trials variable is n. The dependency of the probability of not being ready for rolling is modeled as a function of heating time (Heat), soaking time (Soak), and their interaction. The association parameter displays ordinal measures of association between the observed responses and the predicted probabilities. The ctable parameter saves statistics to the Roc data table that are used to evaluate the predictive power of the model. The output parameter stores the linear predictors and the predicted probabilities in the Out data table along with the other variables from the input data table.

m <- cas.builtins.loadActionSet(s, actionset='regression')

cas.regression.logistic(s,
    association='true',
    ctable=list(casOut=list(name='Roc', replace='true'),
                fpf='FPF',
                nocounts='true',
                tpf='TPF'),
    model=list(depVars=list(list(name='r')),
               effects=list(list(vars=c('Heat', 'Soak')),
                            list(interaction='CROSS', vars=c('Heat', 'Soak'))),
               trial='n'),
    output=list(casOut=list(name='Out', replace='true'),
                copyVars=c('Heat', 'Soak'),
                pred='Pred',
                xBeta='_XBETA_'),
    table='Ingots')

The "Model Information" table shows that the data are modeled as binomially distributed with a logit link function (Output 22.6.34). This is the default link function in the logistic action for binary and binomial data. The action estimates the parameters of the model by a Newton-Raphson algorithm.

The second table in Output 22.6.34 shows that all 19 observations in the input data table were used in the analysis and the "Response Profile" table shows that the total number of events and nonevents equals 12 and 375, respectively.

Output 22.6.35 displays the convergence status table for this run. The logistic action satisfies the gconv convergence criterion.

Output 22.6.36 displays the "Dimensions" table for the model. The design matrix of the model (the matrix) has four columns, which correspond to the intercept, the Heat effect, the Soak effect, and the interaction of the Heat and Soak effects. The model is nonsingular because the rank of the crossproducts matrix equals the number of columns in . All parameters are estimable, and all participate in the optimization.

Output 22.6.37 displays the "Fit Statistics" table for this run. Evaluated at the converged estimates, –2 times the value of the log-likelihood function equals 27.95689. Further fit statistics are also displayed, all of them in "smaller is better" form. The AIC, AICC, and SBC criteria are used to compare non-nested models and to penalize the model fit for the number of observations and parameters. You can use the –2 log-likelihood value to compare nested models by way of a likelihood ratio test.

Output 22.6.38 shows the test of the global hypothesis that the effects jointly do not affect the probability of ingot readiness. You can obtain the chi-square test statistic by comparing the –2 log-likelihood value of the model with covariates to the value in the intercept-only model. The test is significant, with a p-value of 0.0082. One or more of the effects in the model have a significant impact on the probability of ingot readiness.

The "Parameter Estimates" table in Output 22.6.39 displays the estimates and standard errors of the model effects.

Output 22.6.40 displays the "Association Statistics" table, which is produced when you specify the association parameter. The table contains four measures of association for assessing the predictive ability of a model. For more information, see the section Association Statistics.

You can construct the prediction equation of the model from the "Parameter Estimates" table. For example, an observation with Heat equal to 14 and Soak equal to 1.7 has the linear predictor

The probability that an ingot with these characteristics is not ready for rolling is

The logistic action computes these linear predictors and probabilities and stores them in the Out data table.

The following fetch action from the tables action set enables you to display Output 22.6.41, which shows the probability that an ingot with Heat equal to 14 and Soak equal to 1.7 is not ready for rolling:

cas.table.fetch(s,
    format='true',
    index='false',
    table=list(name='Out', where='14 = Heat and 1.7 = Soak'))

The ctable parameter computes statistics for binary response models (on the basis of classifying observations according to whether their predicted probabilities exceed certain values) and stores the results in the Roc data table. For more information, see the section Classification Table and ROC Curves.

The Roc data table contains information you can use to produce an ROC curve like the following:

If you have counts of events and nonevents, then instead of using events/trials syntax, you can use the multinomial-trial syntax to perform the same analysis. In the mycas.Ingots data table, the a variable contains the number of ingots ready for rolling (the nonevents). In the following program, you specify the event variable, followed by the nonevent variable:

cas.regression.logistic(s,
    model=list(depVars=list(list(name='r'), list(name='a')),
               effects=list(list(vars=c('Heat', 'Soak')),
                            list(interaction='CROSS', vars=c('Heat', 'Soak')))),
    table='Ingots')

The only difference from the preceding results is that in this case the "Model Information" table displays the event and nonevent variables; this table is not displayed here.

Binomial data are a form of grouped binary data in which "successes" in the underlying Bernoulli trials are totaled. You can thus unwind data for which you use the events/trials syntax and fit these data by using techniques for binary data.

The following DATA step expands the Ingots data table, which contains 12 events in 387 trials, into a data table that has 387 observations:

cas.dataStep.runCode(s,
    code='
   data Ingots_binary;
      set Ingots;
      do i=1 to n;
        if i <= r then y=1; else y = 0;
        output;
      end;
   run;
   ',
    single='no')

The following statements fit the model that contains the Heat effect, the Soak effect, and their interaction to the data table. The event=’1’ response-variable subparameter in the model parameter ensures that the logistic action models the probability that the variable y takes the value 1.

cas.regression.logistic(s,
    model=list(depVars=list(list(name='y', options=list(event='1'))),
               effects=list(list(vars=c('Heat', 'Soak')),
                            list(interaction='CROSS', vars=c('Heat', 'Soak')))),
    table=list(name='Ingots_binary'))

Output 22.6.43 displays the "Model Information," "Number of Observations," and "Response Profile" tables. The data are now modeled as binary (Bernoulli distributed) with a logit link function. The "Response Profile" table shows that the binary response breaks down into 375 observations where y=0 and 12 observations where y=1.

Output 22.6.44 displays the result of the test of the global null hypothesis and the parameter estimates. These results match those in Output 22.6.38 and Output 22.6.39.

(View the complete code for this example.)

This section uses the data and selected model from Example 22.5 to demonstrate how to request odds ratios. In the following statements, the oddsratio parameter is specified to compute the odds ratios, and the variables for which odds ratios are computed are listed in the vars subparameter.

proc cas;
regression.logistic /
    display="oddsratios",
    model={depvar={{name="y", options={event="1"}}}, effects={"x2", "x8"}},
    oddsratios={vars={"x2", "x8"}},
    table="getStarted";
run;

The resulting odds ratios are displayed in Output 22.7.1.

From Output 22.7.1, subjects with an X2 value that is one unit larger than its current value have 3.163 times the odds of having a Y = 1 response, and subjects with an X8 value that is one unit larger have 2.804 times the odds.

If your model has interactions, then the odds ratio computations depend on the values of the interacting variables. The following code produces odds ratios for the main effects, but because of interactions, these odds ratios depend on the values of the interacting variables as displayed in Output 22.7.2:

proc cas;
regression.logistic /
    display="oddsratios",
    model={depvar={{name="y", options={event="1"}}},
           effects={{interaction="BAR", vars={"x2", "x8"}}}},
    oddsratios={vars={"x2", "x8"}},
    table="getStarted";
run;

The odds ratios in Output 22.7.2 are computed at the mean values of the interacting variables. You can change these values by specifying an at parameter for each individual odds ratio, as follows:

proc cas;
regression.logistic /
    display="oddsratios",
    model={depvar={{name="y", options={event="1"}}},
           effects={{interaction="BAR", vars={"x2", "x8"}}}},
    oddsratios={vars={{at={{value="2", var="x8"}}, var="x2"},
                      {at={{value={5, 10}, var="x2"}}, var="x8"}}},
    table="getStarted";
run;

The following commands customize these odds ratios by using individual options to compute the change in the odds when X2 is increased by 5 and 10 units, and by considering only levels of C equal to "A" or "B" and setting X8 to 1. The results are displayed in Output 22.7.4.

proc cas;
regression.logistic /
    class={"C"},
    display="oddsratios",
    model={depvar={{name="y", options={event="1"}}},
           effects={{interaction="BAR", vars={"C", "x2", "x8"}}}},
    oddsratios={vars={{at={{level={"A", "B"}, var="C"}, {values={1}, var="x8"}},
                       unit={5, 10},
                       var="x2"}}},
    table="getStarted";
run;

You can perform the same customization by specifying global options as follows:

proc cas;
regression.logistic /
    class={"C"},
    display="oddsratios",
    model={depvar={{name="y", options={event="1"}}},
           effects={{interaction="BAR", vars={"C", "x2", "x8"}}}},
    oddsratios={at={{level={"A", "B"}, var="C"}, {values={1}, var="x8"}},
                unit={{value={5, 10}, var="x2"}},
                vars="x2"},
    table="getStarted";
run;

(View the complete code for this example.)

This section uses the data and selected model from Example 22.5 to demonstrate how to request odds ratios. In the following statements, the oddsratio parameter is specified to compute the odds ratios, and the variables for which odds ratios are computed are listed in the vars subparameter.

s:regression_logistic{
   table='getStarted',
   model={depvar={{name='y', options={event='1'}}},
          effects={'x2', 'x8'}},
   oddsratios={vars={'x2', 'x8'}},
   display='oddsratios' }

The resulting odds ratios are displayed in Output 22.7.5.

From Output 22.7.5, subjects with an X2 value that is one unit larger than its current value have 3.163 times the odds of having a Y = 1 response, and subjects with an X8 value that is one unit larger have 2.804 times the odds.

If your model has interactions, then the odds ratio computations depend on the values of the interacting variables. The following code produces odds ratios for the main effects, but because of interactions, these odds ratios depend on the values of the interacting variables as displayed in Output 22.7.6:

s:regression_logistic{
   table='getStarted',
   model={depvar={{name='y', options={event='1'}}},
          effects={{vars={'x2', 'x8'},interaction='BAR'}}},
   oddsratios={vars={'x2','x8'}},
   display='oddsratios' }

The odds ratios in Output 22.7.6 are computed at the mean values of the interacting variables. You can change these values by specifying an at parameter for each individual odds ratio, as follows:

s:regression_logistic{
   table='getStarted',
   model={depvar={{name='y', options={event='1'}}},
          effects={{vars={'x2', 'x8'},interaction='BAR'}}},
   oddsratios={vars={{var='x2',at={{var='x8',value=2}}},
                     {var='x8',at={{var='x2',value={5,10}}}}}},
   display='oddsratios' }

The following commands customize these odds ratios by using individual options to compute the change in the odds when X2 is increased by 5 and 10 units, and by considering only levels of C equal to "A" or "B" and setting X8 to 1. The results are displayed in Output 22.7.8.

s:regression_logistic{
   table='getStarted',
   class={'C'},
   model={depvar={{name='y', options={event='1'}}},
          effects={{vars={'C', 'x2', 'x8'}, interaction='BAR'}}},
   oddsratios={vars={
      {var='x2',unit={5,10},at={{var='C',level={'A','B'}},
                                {var='x8',values={1}} } } }},
   display='oddsratios' }

You can perform the same customization by specifying global options as follows:

s:regression_logistic{
   table='getStarted',
   class={'C'},
   model={depvar={{name='y', options={event='1'}}},
          effects={{vars={'C', 'x2', 'x8'},interaction='BAR'}}},
   oddsratios={vars='x2',
               unit={{var='x2',value={5,10}}},
               at={{var='C',level={'A','B'}},
                   {var='x8',values={1}} }},
   display='oddsratios' }