# How to make a model for the given data in sas?

I have a nonlinear model: $$r10^{\beta_0-\beta_1A-\beta_2B-\beta_4D-\beta_{2,4}BD}$$. I used sas to test different models. Here is the code:

proc nlin data=bcoesdata;
parameters r=10 j=5 k=5 l=1 m=1 p=1;
model y = r*(10**(p-(j*a+k*b+l*d+m*b_d)));
run;


After evaluating simpler models (linear regression did not fit), this one gave the following result: Is the column Approx Pr > F indicating that the model is good since its value is lower than $0.01$? I never worked with nonlinear models before so I do not know what metric is appropriate and I guess that the right metric (test) for goodness of fit probably depends on the function used as a model . However, if there is a way to get a linear model that I have missed, feel free to answer that.

P.S. Data aren't normal, if that counts.

A B C D pid Y


where A, B, C, and D denote four parameters, pid is the id of a unit, and Y is the response variable. Data correspond to a full $2^k$ design. Low and high levels of the parameters are coded with -1 and 1. Each combination of A, B, C, and D values was evaluated on every unit in the sample.

Since the data are not normal, and I could not transform it to normality, I used Friedman's test to determine if changing the level of effect makes significant difference in Y. A, B, D and AD proved to be significant.

EDIT 2: The original question was how to determine goodness-of-fit of a nonlinear model from the beginning.

• By "linear regression did not fit", how did it not fit? Was that linear model based on the logs of the dependent variable with $log_{10}(r)+\beta_0-\beta_1 A-\beta_2 B-\beta_{2,4} BD$ ? – JimB Apr 8 '16 at 17:01
• I got it with proc reg. Its R-square was 0.0169 , adjusted R-square was 0.0031 and Pr > F column was 0.3026. proc reg should inspect the model of the form $\beta_0+\beta_1A+\beta_2B+\beta_{2,4}BD$. – Milos Apr 8 '16 at 17:14
• Yes, but did you take the log of the dependent variable? If you could share your data, that would help pinpoint the issue. You could also construct an R-square from the table you present and that R-square doesn't look very large. – JimB Apr 8 '16 at 17:18
• If I transform the response variable, $y$, to $z=log_{10}{y}$, proc reg returns $0.2534$ for R-squared, $0.2428$ for adjusted R-squared and $<0.0001$ for Pr>F column. So, I inverted this and came up with $10^{...}$ model. I did not expect R-square to be near $1$, so lower values should not represent a big problem. – Milos Apr 8 '16 at 17:21
• I can share the raw data, but later, when I get home. :) – Milos Apr 8 '16 at 17:23

Thanks for sharing the data. Because this is a $2^k$ factorial and there appear to be 18 replications (identified by pid), I would treat this as a mixed model and use PROC MIXED (as you're using SAS) using the log of the y values:

ods graphics on;
proc mixed data=test covtest;
class a b d pid;
model logy = a|b|c|d / solution residual;
random pid;
repeated / group=d;
run;
ods graphics off;


This analyzes the experiment as it was laid out with the addition of estimating separate residual variances for each level of d. (A more complicated variance-covariance structure is likely warranted which I'll touch on later.)

This gives a fixed effect for every combination of the 16 combinations of a, b, c, and d. (I suspect that you mentioned only a, b, d, and b*d because a previous analysis suggested that these were the only significant effects. I'm pretty conservative about these things and would argue that you need to do the analysis as designed rather than performing any variable selection on such designed experiments.)

Here are the tests of each of the 16 combinations: and estimates of variance components and a set of residual plots: In the plot of the studentized residuals vs. the predicted values, there are neighboring pairs of predicted values that seem to have matching residuals which strongly suggests that there's more to the design (either a timing or spatial relationship) than just the sharing of a pid. (This is why some of the "columns" of residual points look a little fuzzy as there are two sets of points that are just slightly offset from each other.) If there is more to the experimental design than just pid, then I'd modify the PROC MIXED statements to account for that.

Estimates of the fixed effects can be obtained using the LSMEANS statement.

Update: I really didn't answer your specific question because I don't see that a nonlinear model is necessary. As a reminder when one is looking for normality, it's not about the dependent variable being approximately normal but rather the residuals. It seems that a linear model using the log of the dependent variable gives an adequate description of your data (other than maybe what appears to be the additional covariance structure). One can obtain the LSMEANS and transform back to the original units. • Thank you for your answer. I have just one more question. Why did you use d in the repeated statement? – Milos Apr 9 '16 at 17:27
• By the way, you are right about the pids. There are three groups of units and I will include this information in the analysis. – Milos Apr 9 '16 at 17:38
• I used group = d because in looking at the simple standard deviations of the 16 combinations, the biggest separation in standard deviations was associated with the two levels of d. The option "group = d" results in two different error variances. This doesn't usually change the cell means much at all but can make a big difference in standard errors of effects and/or tests of hypotheses. – JimB Apr 9 '16 at 23:41