Since it seems to be soapbox time: But saying that it is an industry standard is (to me) like saying there is a tool (say a screwdriver) that I can use to drive screws, open paint cans, chip rust, and turn around to use as a hammer, then I should only learn how to use a screwdriver in new ways, rather than noticing that on the toolbench, that I paid a lot of money for, are a paint lid key, a chisel and a true hammer.  I already own them--why not use them? If I were going to do everything in SQL, I'd do it in something relatively native to the database at hand.  At least, that is the thought of our local IT group. I like Chuck's answer--it is a matter of education, and frankly, sometimes there just isn't time in a production environment to start learning something completely new.  I understand that. Steve Denham ... View more

Just a guess.  With dist=beta, the default link is logit.  You may have structural zeroes.  However, since you have run into the same problem after removing empty cells, I wonder if you are at the other end of the spectrum where ALL of the observations are 1, so that the logit goes off scale.  As I said, this is just a guess.  Depending on the size of the R matrix (number of distinct values of daystoend), you may have some luck if you change the random statement to: random daystoend /residual type=chol subject=id; I fear that the matrix may exceed your memory available though. You might also look at Example 40.4: Quasi-likelihod Estimation for Propotions with Unknown Distribution in the documentation, as another way to attack pseudo-binomial data.  Some of the observations there are zeroes. Steve Denham Message was edited by: Steve Denham ... View more

A comment: Lately there seems to be a lot of questions regarding using PROC SQL to do things that there are other PROCS that are specifically written to accomplish (PROC FREQ, PROC MEANS/SUMMARY, PROC REPORT, etc.), or would be better done inside a DATA step.  Maybe I am too old school, but I just wonder why PROC SQL seems to be the default. Steve Denham ... View more

You may run into difficulty scoring, as I think you might have some problems with the model you are using.  Are x1 and x2 realizations from class1 and class2?  The random statement as it stands fits a random slope model, with a fixed intercept.  Is that correct?  If so, then I suggest the following: proc glimmix data=xxx;       class class1 class2;;        model y(event='1') = time class1 class2 / s dist=binary;       random  time/ subject=pat; run; Now, as to scoring a new dataset... If you had SAS 9.3, I would recommend using the STORE command and PROC PLM.  Since this is 9.2, the best I can suggest is to create the dataset yyy as: data yyy; set xxx dataset_to_be_scored; run; where "dataset_to_be_scored" is the new dataset, but with y set to missing for all observations.  The following should give a dataset with predicted values: proc glimmix data=yyy;       class class1 class2;;        model y(event='1') = time class1 class2 / s dist=binary;       random  time/ subject=pat;      output out=scored pred(blup ilink)=predicted; run; Hope this helps. Steve Denham ... View more

If a sample standard error is zero, by implication the population standard error is zero (why? because it is the only inference you can make about the population parameter), and thus there is no variability in the population.  So any sample will give a CI of zero width. This is one of the killers of small sample size, and why, in our shop, there was a caveat that there would be no analysis unless N>2 for every group.  (Note that this fails to account for internal replication when there are multiple groups.) Steve Denham ... View more

I get really scared (and not by bats) when repeated measures are taken at uneven intervals on subjects like this--things are going to get difficult.  Correct me if I am wrong in my assumption that daystoend is a variable that for any given id specifies something about the timing of when the measurement was taken.  If so, you might consider the following: proc glimmix; class group id daystoend; /* create a continuous version of daystoend that will be used to specify the covariance structure */ dte=daystoend; model behavior1 = group daystoend group*daystoend/ dist=beta; random daystoend /residual type=sp(pow)(dte) subject=id; lsmeans group / diff cl; run; In this code, daystoend may differ for different individuals, depending on the number of sessions an individual is observed.  If I am incorrect in this assumption, we need to find some other variable to index the sessions. Steve Denham ... View more

Some more info about the X and Y variables would be a lot of help.  Looking at the plot that accompanied your post, I didn't see much but noise.  The X variable looks to get very large, so a plot of Y vs. log(X) might be revealing.  I don't think there is going to be anything, but you might look at PROC TRANSREG for Box-Cox transformations as a step towards something. Steve Denham ... View more

On SAS-L, I replied, and since folks don't always read both, I tried to recover it from the archive: This looks like an opportunity to use PROC GLIMMIX, and use of the LINK option. For depth LINK=LOG, and for mass LINK=POWER(0.5). Then in the LSMEANS statement, use the ILINK option, and the final values will include the estimates and their standard errors on both the transformed and original scale. The documentation reports that the standard errors on the inverse lined scale are computed by the delta method. I hope this helps. Steve Denham Message was edited by: Steve Denham ... View more

Or, since it looks like there is some variability accounted for by intercept study (0.000429), you might try: proc glimmix data=count_data; class study; model tox =  / dist=nb link=log offset=logPA2 solution; random intercept/subject=study: random intercept/subject=study group=study type=chol(1) solution; parms   (0)            (0.125) (0.333) (0.031) (0.2) (0.0625)                               (0.0169) (0.005) (0.007) (0.055) (0.0476)            (0.083) (0.033) (0.1)  (0.0625)  (0.0212)           (1) /  hold=(2 to 16);  run; This will estimate the first and last parameter, and hold the parameters for the fifteen studies constant.  I would still stay with the Cholesky parameterization of the G matrix. Steve Denham ... View more

I worry about the three parms instead of 20 error. I think it is coming from the following: there are three parms: first random statement = variance between studies, second random statement = also variance between studies, and residual/scale. This would also explain why you get G matrix not positive definite, as the matrix contains a redundant column when both random statements are included. What about: proc glimmix data=count_data; class study; model tox =  / dist=nb link=log offset=logPA2 solution; random intercept/subject=study group=study type=chol(1) solution; parms              (0.125) (0.333) (0.031) (0.2) (0.0625)                               (0.0169) (0.005) (0.007) (0.055) (0.0476)            (0.083) (0.033) (0.1)  (0.0625)  (0.0212)           (1) /  hold=(1to 15);  run; Going to a single random statement and using the Cholesky root.  I offer type=chol, as it ensures that the varainace-covariance matrix is at least positive semidefinite, so that should help, as will removing the redundant random intercept/subject=study line. Good luck. Steve Denham ... View more

I've never attempted anything like this, but I did notice one curiosity in your code for the second model.  You have two random statements, but they are much the same, so that the first statement results in an odd occurrence.  What about a statement that looks like: random intercept/subject=study type=un(1): I think this would give a diagonal matrix, with the diagonal having the variances for each study.  Perhaps lvm will comment, and help out on this one. Steve Denham ... View more