Problem: We have observed erratic behavior of a split-plot model when run with different subsets of our data. For two of six response-variable columns, the script produces incorrect outputs for the random factor in the model. Specifically, the random factor accounts for 0% of the total variance in the REML Variance Component Estimates table, and the table of Least Squares Means displays a column of means but not of standard errors.
We have checked the column properties and replaced the data in the problematic columns, to no avail. Two lines of evidence suggest that this behavior is not due to a script error:
- The script works properly on the other four of response variable columns
- A reduced model—including only the random factor and the problematic data columns—produces incorrect output for one of the response variables while for the other variable the output appears to be correct
We include a brief summary of our experimental design and approach to analysis, as well as examples of the problematic data, the scripts used, and model output. We also attach a pdf showing the script and output for the reduced models described above.
Thanks for any assistance in solving this puzzle!
Experimental design: An ecological field experiment consisting of four (or three) spatially defined blocks (replicates) to which 3 factors were randomly assigned in a split-plot design. At the whole plot scale,
two two-level factors—Limpet removal (+/-) and Moisture addition (+/-), both relatively hard to change—were randomly assigned to four plots within each block. At the sub-plot scale, Blade removal (Blades +/-) was randomly assigned to two sub-plots within each whole plot (see sample layout, below).
Sample dataset
Blk# | Limp | Moist | Blades | Ja88Ap88 | Mr89Nv89 |
1 | +L | +M | +B | 0 | -15 |
1 | -L | +M | +B | -5.88235 | -2.94118 |
1 | +L | -M | +B | 5.88235 | -11.7647 |
1 | -L | -M | +B | -14.4118 | -38.2353 |
2 | +L | +M | +B | -2.94118 | 2.94118 |
2 | -L | +M | +B | 0 | -20.5882 |
2 | +L | -M | +B | 0 | -15 |
2 | -L | -M | +B | 0.29412 | -2.94118 |
3 | +L | +M | +B | 0 | . |
3 | -L | +M | +B | -32.3529 | . |
3 | +L | -M | +B | -17.6471 | . |
3 | -L | -M | +B | -15 | . |
4 | +L | +M | +B | 0 | -2.64706 |
4 | -L | +M | +B | -5.88235 | -14.4118 |
4 | +L | -M | +B | -5.58824 | -43.8235 |
4 | -L | -M | +B | -5.88235 | -35.2941 |
1 | +L | +M | -B | -5.88235 | 11.47059 |
1 | -L | +M | -B | 3.23529 | 5.29412 |
1 | +L | -M | -B | -2.94118 | -8.82353 |
1 | -L | -M | -B | -2.94118 | -11.7647 |
2 | +L | +M | -B | -5.88235 | 6.17647 |
2 | -L | +M | -B | 0 | -6.17647 |
2 | +L | -M | -B | 2.94118 | -26.4706 |
2 | -L | -M | -B | 3.23529 | -20.2941 |
3 | +L | +M | -B | -14.7059 | . |
3 | -L | +M | -B | 14.70588 | . |
3 | +L | -M | -B | -0.29412 | . |
3 | -L | -M | -B | -5.88235 | . |
4 | +L | +M | -B | -8.82353 | -26.4706 |
4 | -L | +M | -B | -26.4706 | -20.5882 |
4 | +L | -M | -B | -14.4118 | -44.1177 |
4 | -L | -M | -B | 2.94118 | -47.3529 |
Analysis. To simplify the analysis, we chose to model the four whole-plot treatment combinations as 4 levels of a single factor (Concatenate[Limp,Moist] in the script, below). We plan to use post-hoc tests to distinguish among treatment means. The random effect (error term for the whole plot) is a concatenation of the block, limpet, and moisture codes (i.e. :Name( "Concatenate[Blk#,Limp,Moist]" )[:Blk#], below).
Sample script for full model.
</CODE>Fit Model(<BR />Y( :Ja88Ap88 ),<BR />Effects(<BR />:Blk#,<BR />:Name( "Concatenate[Limp,Moist]" ),<BR />:Blades,<BR />:Blades * :Name( "Concatenate[Limp,Moist]" )<BR />),<BR />Random Effects( :Name( "Concatenate[Blk#,Limp,Moist]" )[:Blk#] ),<BR />NoBounds( 1 ),<BR />Personality( "Standard Least Squares" ),<BR />Method( "REML" ),<BR />Emphasis( "Minimal Report" ),<BR />Run(<BR />)<BR />),</PRE><P> </P><P><U>Partial output for full model</U>.</P><P><STRONG>REML Variance component estimates </STRONG>For Ja88Apr88 (n = 4 blocks)</P><TABLE><TBODY><TR><TD><P><STRONG>Random Effect</STRONG></P></TD><TD><P><STRONG>Var Ratio</STRONG></P></TD><TD><P><STRONG>Var Component</STRONG></P></TD><TD><P><STRONG>Std Error</STRONG></P></TD><TD><P><STRONG>95% Lower</STRONG></P></TD><TD><P><STRONG>95% Upper</STRONG></P></TD><TD><P><STRONG>Wald p-Value</STRONG></P></TD><TD><P><STRONG>Pct of Total</STRONG></P></TD><TD> </TD></TR><TR><TD><P>Concatenate[Blk#,Limp,Moist][Blk#]</P></TD><TD><P> -0.324822</P></TD><TD><P> -40.38441</P></TD><TD><P>27.376421</P></TD><TD><P> -94.04121</P></TD><TD><P>13.272392</P></TD><TD><P>0.1402</P></TD><TD><P>0.000</P></TD><TD> </TD></TR><TR><TD><P>Residual</P></TD><TD><P> </P></TD><TD><P>124.32778</P></TD><TD><P>50.756605</P></TD><TD><P>63.930876</P></TD><TD><P>338.78407</P></TD><TD><P> </P></TD><TD><P>100.000</P></TD><TD> </TD></TR><TR><TD><P>Total</P></TD><TD><P> </P></TD><TD><P>124.32778</P></TD><TD><P>50.756605</P></TD><TD><P>63.930876</P></TD><TD><P>338.78407</P></TD><TD><P> </P></TD><TD><P>100.000</P></TD><TD> </TD></TR></TBODY></TABLE><P> </P><P><STRONG>Fixed effect tests </STRONG>For Ja88Apr88 <U>Note</U> denominator degrees of freedom appear to be correct.</P><TABLE><TBODY><TR><TD><P><STRONG>Source</STRONG></P></TD><TD><P><STRONG>Nparm</STRONG></P></TD><TD><P><STRONG>DF</STRONG></P></TD><TD><P><STRONG>DFDen</STRONG></P></TD><TD><P><STRONG>F Ratio</STRONG></P></TD><TD><P><STRONG>Prob > F</STRONG></P></TD><TD> </TD></TR><TR><TD><P>Blk#</P></TD><TD><P>3</P></TD><TD><P>3</P></TD><TD><P>9</P></TD><TD><P>3.1201</P></TD><TD><P>0.0808</P></TD><TD> </TD></TR><TR><TD><P>Concatenate[Limp,Moist]</P></TD><TD><P>3</P></TD><TD><P>3</P></TD><TD><P>9</P></TD><TD><P>0.2215</P></TD><TD><P>0.8791</P></TD><TD> </TD></TR><TR><TD><P>Blades</P></TD><TD><P>1</P></TD><TD><P>1</P></TD><TD><P>12</P></TD><TD><P>0.3675</P></TD><TD><P>0.5557</P></TD><TD> </TD></TR><TR><TD><P>Blades*Concatenate[Limp,Moist]</P></TD><TD><P>3</P></TD><TD><P>3</P></TD><TD><P>12</P></TD><TD><P>1.0059</P></TD><TD><P>0.4238</P></TD><TD> </TD></TR></TBODY></TABLE><P> </P><P><STRONG>Least Squares Means Table</STRONG> For Ja88Apr88 (n = 4 blocks)</P><TABLE><TBODY><TR><TD><P><STRONG>Level</STRONG></P></TD><TD><P><STRONG>Least Sq Mean</STRONG></P></TD><TD><P><STRONG>Std Error</STRONG></P></TD><TD> </TD></TR><TR><TD><P>[1]1, +L, +M</P></TD><TD><P> -2.28820</P></TD><TD><P>.</P></TD><TD> </TD></TR><TR><TD><P>[1]1, +L, -M</P></TD><TD><P> -9.03709</P></TD><TD><P>.</P></TD><TD> </TD></TR><TR><TD><P>[1]1, -L, +M</P></TD><TD><P> -8.62806</P></TD><TD><P>.</P></TD><TD> </TD></TR><TR><TD><P>[1]1, -L, -M</P></TD><TD><P>8.48276</P></TD><TD><P>.</P></TD><TD> </TD></TR><TR><TD><P>[2]2, +L, +M</P></TD><TD><P>7.78410</P></TD><TD><P>.</P></TD><TD> </TD></TR><TR><TD><P>[2]2, +L, -M</P></TD><TD><P> -1.69162</P></TD><TD><P>.</P></TD><TD> </TD></TR><TR><TD><P>[2]2, -L, +M</P></TD><TD><P> -3.73673</P></TD><TD><P>.</P></TD><TD> </TD></TR><TR><TD><P>[2]2, -L, -M</P></TD><TD><P> -3.53222</P></TD><TD><P>.</P></TD><TD> </TD></TR><TR><TD><P>[3]3, +L, +M</P></TD><TD><P> -11.31712</P></TD><TD><P>.</P></TD><TD> </TD></TR><TR><TD><P>[3]3, +L, -M</P></TD><TD><P> -6.88603</P></TD><TD><P>.</P></TD><TD> </TD></TR><TR><TD><P>[3]3, -L, +M</P></TD><TD><P> -11.93065</P></TD><TD><P>.</P></TD><TD> </TD></TR><TR><TD><P>[3]3, -L, -M</P></TD><TD><P> -5.45444</P></TD><TD><P>.</P></TD><TD> </TD></TR><TR><TD><P>[4]4, +L, +M</P></TD><TD><P> -14.25232</P></TD><TD><P>.</P></TD><TD> </TD></TR><TR><TD><P>[4]4, +L, -M</P></TD><TD><P> -2.45880</P></TD><TD><P>.</P></TD><TD> </TD></TR><TR><TD><P>[4]4, -L, +M</P></TD><TD><P>4.22192</P></TD><TD><P>.</P></TD><TD> </TD></TR><TR><TD><P>[4]4, -L, -M</P></TD><TD><P> -19.56963</P></TD><TD><P>.</P></TD><TD> </TD></TR></TBODY></TABLE><P> </P></LINGO-BODY><LINGO-LABS id="lingo-labs-328571" slang="en-US"><LINGO-LABEL>Data Table</LINGO-LABEL><LINGO-LABEL>Scripting</LINGO-LABEL></LINGO-LABS><LINGO-SUB id="lingo-sub-328695" slang="en-US">Re: Erratic behavior of split-plot model: Contaminated data, script error, or bug?</LINGO-SUB><LINGO-BODY id="lingo-body-328695" slang="en-US"><P>Let me preface my comments by I don't understand the situation you are describing from a scientific or engineering point of view so my comments may be completely off base. It does make sense that changing limpets and moisture would be difficult to do randomly. I don't know what your hypotheses are nor do I understand the response variable or the factors for that matter? I have some questions regarding your experiment. The difficulty I am having is understanding your design plan (based on the sample data you provided):</P><P>1. This statement makes no sense to me: "<SPAN>to which 3 factors were randomly assigned in a split-plot design". Split-plot designs purposely restrict randomization, so I don't understand what you are doing. In the sample data set, it looks like your "hard to change f</SPAN>actors" are changing a lot! The subplot factor doesn't appear to change nearly as much as your whole plot factors? Do you have the run order?</P><P>2. I also don't understand your use of the block? I presume you are using the blocks to get an estimate of the whole plot errors? Let me see if I understand what you are doing...You are running a 2^3 factorial experiment in 4 different "locations/spaces" (blocks). Within each block you are intending to restrict randomization on 2 of 3 factors (L and M) in a factorial matrix. This sounds more like a split-split-plot. I added just the design in a jmp file to see if this is what your are running.</P><P>3. Why would you combine L & M? If you do this you lose the ability to quantify the L*M interaction and any other 2nd order effect with L and M? Instead you can estimate quadratic and cubic effects??? </P><P>4. Your measurement system look to have a discrimination issue (it doesn't appear continuous). There are only some numbers available looking at the sample data set? Is this the real data?</P><P>Food for thought (pardon my rambling)... the reason you need to be aware and take into account restrictions on randomization is due to the potential for unanticipated changes in the noise over the course of the experiment which in turn impacts the experimental errors and subsequently the statistical tests. The thinking is the whole plot errors will be different than the sub-plot errors and therefore a different test of significance should be performed for each. If in planning your experiment, if you identify the noise variables and hypothesize about the affects of those noise variables, you can take that into account in the designing of your experiment (e.g., what is changing within block and what is changing between block?). If there is very little difference between the noise changing at the whole plot level (between block) and the noise changing at the sub-pot level, then taking into account the restriction may be of no consequence. If, however, there are significant differences in the noise (and therefore the estimated errors) at the whole plot level and the noise at the sub-plot level, then tests of significance should reflect that. </P><P> </P></LINGO-BODY><LINGO-SUB id="lingo-sub-328867" slang="en-US">Re: Erratic behavior of split-plot model: Contaminated data, script error, or bug?</LINGO-SUB><LINGO-BODY id="lingo-body-328867" slang="en-US"><P>Please read <A href="https://community.jmp.com/t5/user/viewprofilepage/user-id/4358" target="_blank">@statman</A>'s reply carefully.</P>
<P> </P>
<P>I am also confused by your post:</P>
<UL>
<LI>You have concatenated the Limpet and Moisture levels and then treated their joint changes as a random effect on the response? Wouldn't you expect these factors to produce fixed effects?</LI>
<LI>What you call Blk# is actually the whole plot? That is, these runs reset the hard to change factor once? Then it should probably be the only random effect in the model.</LI>
<LI>It appears that Blades is the whole plot (hard to change) factor, but then how is Blades set to both -B and +B in the same whole plot?</LI>
<LI>Your script indicates that you are using the same factors for both the fixed and random effects. Is that intentional?</LI>
<LI>You have lost all of whole plot = 3?</LI>
</UL>
<P>Perhaps this model will provide the information that you want as more stable estimates:</P>
<P> </P>
<PRE><CODE class=" language-jsl">Fit Model(
Y( :Ja88Ap88 ),
Effects( :Limp, :Moist, :Blades, :Blades * :Limp, :Blades * :Moist ),
Random Effects( :Blk# ),
No Bounds( 1 ),
Personality( "Standard Least Squares" ),
Method( "REML" ),
Emphasis( "Minimal Report" ),
Run()
);
