# How to analyze this messy design?

I need to analyze a data-set, with a very messy design, I am not sure how. I will try to make it simple. A new kind of stitches was invented, and is tested vs. 2 old kind of stitches. I will call this: Treatment, control 1 and control 2.

A few patients were selected, needing stitches. Each patient was assigned (by the need, not by random) to one of two stitching techniques. There are others, but only these two were tested. I call this variable "procedure type".

Each patient, needed one stitch or more, up to 5 stitches. The problem is, that those patients in need for more than one stitch, received (sometimes) a mix of different stitches. Meaning, that one patient could have treatment and control 1, and other could have 2 treatment and 2 control 2, etc...

I hope I am not mixing terms here, but my experimental unit is the patient and the observations unit is the stitch within the patient. I have 40 treatment stitches, and 20 of each control, all coming from 25 patients.

The measure being tested, is how many minutes did it take for each stitch until the bleeding stopped. The variable can take values: 0,0.5,1,1.5,2,2.5,...

Clinicians claim that the procedure type, i.e. the technique being used is of no clinical importance and should not have any affect. Thus they think the analysis could ignore this factor.

I wanted to ask you, how would you analyze this data, with and without taking into account the technique. I think it should be some sort of a mixed model or a generalized estimating equations model, but I am not sure to define this design exactly in terms of "nesting", "blocking" and so on, and I would like to also fit a SAS code to model the differences between the means of the 3 different treatments, and I wouldn't mind verifying it with R.

For your convenience, I attach a diagram I made of the design: the blue circles are treatment stitches, the red ones are control 1 stitches, and the green are control 2 stitches. Each yellow rectangle is a patient, and the big brown rectangles are the techniques, which again, I would like to try to fit and to try to ignore the see the outcome of that.

• Am I correct that the closest known design is the split plot ? Aug 17, 2014 at 10:36
• It could be a split plot, but as you do not care about procedure type (and that it is not suppose to change anything), I think you could mix them together and take them out of your design... Aug 18, 2014 at 12:59

Since the procedure type is of no clinical importance, it can be ignored or included as a random factor to see if it is statistically significant. So the random factors would be either patients only or patients nested within the procedures.

The dependent variable can be modeled as the geometric distribution because it is the number of minutes it took until the first success in stopping bleeding. You can specify it as a negative binomial distribution with $\theta = 1$. I know that this is an unbalanced design, the number and types of stitches vary across patients, but (restricted) maximum likelihood estimation should be able to estimate population parameters accurately, assuming that stitch types and numbers are distributed across the patients and procedures (completely) randomly.

Sample size permitting, here are some other random effects you may want to include in the model are: the order in which stitching was performed within each patient and surgeons who performed stitches (if there were two or more).

I am sure that you can model mixed effects negative-binomial regression with SAS, and I know you can with R (e.g., the glmmADMB or lme4 package).

• Since the DV takes non integer values in the sample data, is it really wise to model it as geometric? (I understand the argument for "time until first success", though it isn't clear to me if there would be independence between time periods, but the time interval used was 30 seconds so perhaps the data should be doubled first?) Nov 7, 2014 at 11:48

If you have enough patients, perhaps you could ignore the first set of stitches, and analyze only the time to stop bleeding after the second set. Then you might categorize the patients into four groups: received all control stitches, received all control 1 stitches, received all control 2 stitches, received all treatment stitches, received some mixture of stitches. Then see if there is any diference in mean time between these four groups. If so, then look at your data further, perhaps looking at other variables, including the first set of stitches. Also, perhaps the time to stop bleeding after the first set of stitches could be analyzed with the same categorization. BTW, there seems to be more to the situation than you describe, or why would a clinician use 2 different types of stitches in either set?

Suppose your data are in a SAS data set called Stitches and are organized like this:

Procedure  Patient     Stitch  Time
1       1  Treatment   0.0
1       1  Treatment   1.0
1       1  Control-1   1.5
1       1  Control-1   0.0
1       2  Treatment   1.0
1       2  Control-1   0.0
1       3  Treatment   1.5
1       3  Control-2   3.0
...
2      25  Treatment   0.0
2      25  Control-1   1.0
2      25  Control-1   0.5


Then, here is some untested SAS code that will fit a very minimal linear model:

proc mixed data=Stitches;
class Procedure Stitch Patient;
model Time = Procedure Stitch Procedure*Stitch;
random Patient(Procedure);
run;


If you have a lot of zeroes, this probably isn't going to fit well. If the data distribution is too granular, this probably isn't going to fit well. It's probably not going to fit well, is what I'm saying.

I definitely agree with Masato Nakazawa that it would be nice to incorporate a random effect for surgeon if that makes sense, and also a fixed effect representing stitch order. You could start with stitch order as a simple covariate (i.e., don't put it in the CLASS statement).

Again, untested code, but you can morph the code above in SAS to use a generalized linear mixed effect model and either use Masato Nakazawa's suggestion of negative binomial link or a Poisson link as a starting point:

proc glimmix data=Stitches;
class Procedure Stitch Patient;
model Time = Procedure Stitch Procedure*Stitch / dist=poisson link=log;
random Patient(Procedure);
run;


If appropriate, you could use least squares means as a first cut to evaluating differences among stitch treatments:

  lsmeans Stitch / cl pdiff ilink;


Or, you could set up a contrast between the new treatment and the average of the two controls.

I would keep the procedure type in the model. After all, if it is truly unimportant, then it isn't going to hurt much to lose a few degrees of freedom.

There are many issues that may crop up for you, though. With the linear model approach, there's most likely going to be some fit problems. It might still be pragmatic to use that approach. With the generalized linear model approach, you may have fitting problems.

If the interaction of procedure and stitch is significant, then you should probably slice out your results by procedure.

For model structure, you could make arguments that the stitches of the same type within a patient might be better considered sub-samples. Or, you might have position effects as well as stitch order effects. Right now, this model just assumes that patients differ in their clotting time.

Model diagnostics should be assessed. The data should be plotted.

Again, all this code is untested, as I don't have SAS on this machine right now! I hope someone will correct me if I mis-specified it!