procD.lm and RRPP - how to disable model truncating

Tyler Serdar

unread,

Feb 23, 2022, 4:51:36 PM2/23/22

to geomorph R package

Dear all,

I have set up a fairly complex procD.lm/RRPP model that contains a three-way interaction and a random effect.

Now, when I add the random effect to the model, the automatic model truncating removes the three-way interaction as well as one of the two-way interactions from the "full model" results.

However, when I do not add the random effect to the model, all interactions are included in the full model.

Now, the reviewers are asking me to a) include the random effect in any case to control for the related variation, and b) to report the full model entirely (i.e. including the three-way interaction and all two-way interactions).

Thus, I would like to ask you if you know any way of turning off model truncating in procD.lm/RRPP so I can actually have a complete "full model" to report. Even if some interactions are redundant, I need to be able to report the (likely extremely small) effect size for every single interaction.

Thanks for your help!

Tyler Serdar

Mike Collyer

unread,

Feb 26, 2022, 7:15:16 PM2/26/22

to geomorph R package

Tyler,

This sounds like an R formula issue (procD.lm does not actively truncate model terms). Would you mind sharing your formulae so it is obvious what you mean by adding random effects?

Thanks,

Mike

--
You received this message because you are subscribed to the Google Groups "geomorph R package" group.
To unsubscribe from this group and stop receiving emails from it, send an email to geomorph-r-pack...@googlegroups.com.
To view this discussion on the web, visit https://groups.google.com/d/msgid/geomorph-r-package/5051dab6-a3fe-472d-a3d4-cf2dafd0af09n%40googlegroups.com.

Tyler Serdar

unread,

Feb 28, 2022, 3:14:54 AM2/28/22

to geomorph R package

Dear Mike,

m.1<-procD.lm(coords~a*b*c, iter=9999, data=dat.full2, RRPP=FALSE, SS.type="III", print.progress=FALSE) anova(m.1)

is working.

m.1<-procD.lm(coords~a*b+family, iter=9999, data=dat.full2, RRPP=FALSE, SS.type="III", print.progress=FALSE) anova(m.1)

is working.

However,

m.1<-procD.lm(coords~a*b*c+family, iter=9999, data=dat.full2, RRPP=FALSE, SS.type="III", print.progress=FALSE) anova(m.1)

does not work, and errors with:

ERROR: 'names' attribute [8] must be the same length as the vector [7]

If I run the same model with RRPP or procd.lm and RRPP=TRUE, then it gets truncated:

m.1<-procD.lm(coords~a*b*c+family, iter=9999, data=dat.full2, RRPP=TRUE, SS.type="III", print.progress=FALSE) anova(m.1)

Warning: Because variables in the linear model are redundant,
the linear model design has been truncated (via QR decomposition).
Original X columns: 50
Final X columns (rank): 47
Check coefficients or degrees of freedom in ANOVA to see changes.

outputted effects: (note: a*b is missing!)

a

b

c

family

b*c

a*c

a*b*c

If I run:

m.1<-procD.lm(coords~a*b*c*family, iter=9999, data=dat.full2, RRPP=TRUE, SS.type="III", print.progress=FALSE) anova(m.1)

WARNING: Because variables in the linear model are redundant,
the linear model design has been truncated (via QR decomposition).
Original X columns: 69
Final X columns (rank): 49
Check coefficients or degrees of freedom in ANOVA to see changes.

outputted effects: (note: a*b and c*family are missing!)

a

b

c

family

a*c

b*c

a*b*c

Mike Collyer

unread,

Feb 28, 2022, 6:04:26 AM2/28/22

to geomorph R package

Dear Tyler,

When you see the warning messages, it’s because the linear model design matrix is not full rank. RRPP is not removing terms from the model but QR decomposition of the linear model design matrix results in fewer columns and it is possible terms are removed by culling the columns. What makes the design matrix lose rank? There are two possible reasons this happens. One reason is that you have more independent variables than observations. Here is an example:

> y <- rnorm(3)

> x1 <- rnorm(3)

> x2 <- rnorm(3)

> x3 <- rnorm(3)

> x4 <- rnorm(3)

>

> X <- model.matrix(y ~ x1 + x2 + x3 + x4)

> X

(Intercept) x1 x2 x3 x4

1 1 0.8556341 -1.1927499 0.03386379 0.5490253

2 1 -0.4750418 -0.5066449 -0.33475282 1.0312301

3 1 -0.5038879 -1.1455421 1.31046716 -0.8673995

attr(,"assign")

[1] 0 1 2 3 4

> dim(X)

[1] 3 5

> qr(X)$rank

[1] 3

The rank of the design matrix is defined by the number of rows rather than the number of columns, so the matrix will be reduced to have only 3 columns. The second reason is that variables are redundant. another example:

> x2 <- x1* 5

> X2 <- model.matrix(y ~ x1 + x2 )

> X2

(Intercept) x1 x2

1 1 0.8556341 4.278171

2 1 -0.4750418 -2.375209

3 1 -0.5038879 -2.519439

attr(,"assign")

[1] 0 1 2

> dim(X2)

[1] 3 3

> qr(X2)$rank

[1] 2

So in this case, a linear model design attempts to have two variables, but since one variable is a linear transformation of the other, the columns are reduced.

It looks like your attempt suffers from the former issue, because anova is expecting 8 terms but finding only 7. There is nothing that can be done to not “truncate" the linear model design matrix. There is no way to create sampling distributions for statistics based on coefficients that do not exist, nor can exist. If we did not use QR decomposition to fix the linear model design matrix, it would not be possible to estimate coefficients because that would involve trying to invert singular matrices. It simply does not work.

You either need a simpler model design or more observations (research subjects) to make it work.

Mike

To view this discussion on the web, visit https://groups.google.com/d/msgid/geomorph-r-package/8c45c167-d03d-43a5-9437-d11c163657d6n%40googlegroups.com.

Tyler Serdar

unread,

Feb 28, 2022, 6:19:44 AM2/28/22

to geomorph R package

Dear Mike,

thank you for your advice.

the issue is likely that the random factor "family" is not fully factorial, so when I add it, not every treatment combination is present within every family. The sample size is around 3000, but individual families can contribute only 10-20 samples. With my previous experience, classical mixed-effect models (lme) do not suffer from this issue, which is why I can add random factors that are not fully factorial. What would you advise for me to do to account for the random factor "family" here, given that family sample sizes are too low?

So far, I only see 2 options:

a) Run the model without adding the family random effect.

b) Run the model with family as the random factor, and be unable to report effect sizes for some interactions.

My main gripe is that both approaches have (vastly) different effect sizes for a, b, and c.

Cheers,

Tyler

Mike Collyer

unread,

Feb 28, 2022, 6:59:20 AM2/28/22

to geomorph R package

Although I am mostly guessing — without your actual data — I think that because a * b is missing in various attempts, it suggests that this interaction is invariant; i.e., certain levels of b are only found within certain levels of a, or vice versa. This might be overcome by (my original hypothesis) changing the formula, such that it has a/b or b/a rather than a*b. This is just a guess.

But assuming it is caused potentially by the family effect, then what is your purpose? To control for family variation? It is possible to perform an lme-like analysis assuming that least squares estimates of coefficients are sufficient by doing the following:

fit.rand <- lm.rrpp(Y ~ family, …)

anova(fit.rand) # evaluate random effect

fit.null <- fit.rand <- lm.rrpp(Y ~ 1, …)

# remove family effect

Yf <- fit.null$LM$fitted + fit.rand$LM$residuals

# run a fixed-effects only model, having held constant the random effects

fit.fixed <- lm(Yf ~ a * b * c, …)

anova(fit.fixed)

If you do this and a * b is still missing, at least you will know why.

Hope that helps,

Mike

To view this discussion on the web, visit https://groups.google.com/d/msgid/geomorph-r-package/98b4c68c-ed92-459c-835c-9e1d12d29c16n%40googlegroups.com.

Tyler Serdar

unread,

Mar 2, 2022, 6:33:42 AM3/2/22

to geomorph R package

Dear Mike,

thank you for your advice. Your proposed way of controlling for family variation indeed allows me to run the full model, and I do get results for all interaction terms, including a*b.

However, given that the model results are vastly different from the ones without this random effect, there is just one more thing I'd like to ask you for clarification. My fixed effects a and b can never share the same families. The same families can only be present within the fixed effect c. I wonder if this fact leads to incorrect model results, and whether you'd suggest me to control for family variation in a different way given this knowledge.

Thanks so much for your help again!

Tyler

Mike Collyer

unread,

Mar 2, 2022, 7:49:55 AM3/2/22

to geomorph R package

Hi Tyler,

If you have families that are unique to factor a and families that are unique to factor b, then factors a and b are confounded with family. There is no easy way to test effects for a and b, independent of family, and it seems a three-way interaction might not be possible.

I’m having trouble understanding, e.g., an experiment, where certain treatments do not interacte. If this is the case, then having and a * b * c interaction in the model makes no sense and would be why a * b disappears. an experiment that only has a * c and b * c interactions would make it possible for each family to be represented in c, yet be unique to a and b. Removing the family effect as I suggested might work for you, analytically, but under the mistaken assumption that you rendered the observations independent.

Maybe you can explain these details better (what are the factor levels in a, b, and c, and how it is possible that certain interactions are possible while others are not?). I’m thinking the reviewers might have sent you down a rabbit hole that does not have the solution for what is requested.

Mike

To view this discussion on the web, visit https://groups.google.com/d/msgid/geomorph-r-package/9854a33e-ff7c-4fba-a8ef-aa7bc1f33470n%40googlegroups.com.

Message has been deleted

Tyler Serdar

unread,

Mar 2, 2022, 11:32:09 AM3/2/22

to geomorph R package

Dear Mike,

okay - let me try to explain it better.

The experiment is based on two distinct environments. These two environments were crossed in mothers, fathers, and offspring in a 2 x 2 x 2 design, resulting in total 8 treatments.

Hence, the a*b*c interaction is actually a maternal environment * paternal environment * offspring environment interaction.

Now, one mother and one father, a single family, generated many offspring. These offspring were then split up into the two environments. Hence, the same family can be present across offspring environments.

However, mothers and fathers could not be re-used given the experimental design. Hence, every "family" exists only once within each maternal-paternal environment combination (4 different combinations).

I hope this makes you understand the challenge better - it certainly feels like a rabbit hole to control for family effects while determining the effects of maternal, paternal and offspring environments and their interactions.

Thank you for your help again!

Tyler

Mike Collyer

unread,

Mar 3, 2022, 7:36:15 AM3/3/22

to geomorph R package

Hi Tyler,

Maybe somebody who is more experienced with quantitative genetics designs can offer a better explanation than me. I believe it is probably a formula problem with you trying to define interactions, when some interactions between mothers and fathers (family levels) do not exist in all environments. If I understand the design correctly, I can simulate it and do a practice analysis like this (note that the error terms for F stats in ANOVA are not adjusted):

> father <- factor(rep(letters[1:4], each = 4))

> mother <- factor(rep(letters[5:8], each = 4))

> A <- factor(rep(1:2, each = 8))

> B <- factor(rep(rep(1:2, each = 4), 2))

>

> C <- factor(rep(1:2, 8))

> y <- rnorm(16)

>

> df <- data.frame(A = A, B = B, father = father, mother = mother, C = C, y = y)

> df

A B father mother C y

1 1 1 a e 1 -0.854084818

2 1 1 a e 2 -0.573175241

3 1 1 a e 1 -1.262719277

4 1 1 a e 2 -0.644844684

5 1 2 b f 1 0.045082777

6 1 2 b f 2 0.829667393

7 1 2 b f 1 -0.005103052

8 1 2 b f 2 0.240001039

9 2 1 c g 1 -0.192067928

10 2 1 c g 2 -0.505950086

11 2 1 c g 1 0.841230889

12 2 1 c g 2 0.187821144

13 2 2 d h 1 -0.454320631

14 2 2 d h 2 -0.760631043

15 2 2 d h 1 0.926958511

16 2 2 d h 2 -0.028978510

>

> fit <- lm.rrpp(y~ A / (father:mother) + B / (father:mother) + A * B * C, data = df)

Warning: Because variables in the linear model are redundant,

the linear model design has been truncated (via QR decomposition).

Original X columns: 56

Final X columns (rank): 8

Check coefficients or degrees of freedom in ANOVA to see changes.

> anova(fit)

Warning ANOVA is missing some terms, likely because

some independent variables were redundant.

If the residual SS is 0, results should not be trusted

Analysis of Variance, using Residual Randomization

Permutation procedure: Randomization of null model residuals

Number of permutations: 1000

Estimation method: Ordinary Least Squares

Sums of Squares and Cross-products: Type I

Effect sizes (Z) based on F distributions

Df SS MS Rsq F Z Pr(>F)

A 1 0.3134 0.31339 0.05061 1.1107 0.53677 0.298

B 1 0.9008 0.90082 0.14549 3.1926 1.22115 0.114

C 1 0.0057 0.00567 0.00091 0.0201 -1.25962 0.886

A:B 1 1.6208 1.62083 0.26178 5.7443 1.71771 0.034 *

A:C 1 1.0806 1.08057 0.17452 3.8296 1.28432 0.095 .

B:C 1 0.0017 0.00168 0.00027 0.0060 -1.66583 0.948

A:B:C 1 0.0113 0.01133 0.00183 0.0402 -1.04955 0.843

Residuals 8 2.2573 0.28216 0.36458

Total 15 6.1916

---

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Call: lm.rrpp(f1 = y ~ A/(father:mother) + B/(father:mother) + A * B * C, data = df)

>

So that works. Notice that nested terms are not part of the ANOVA. Then, if I create a family variable as I think you describe it, where not all families are represented in every level of A or B, I get this:

> family <- interaction(father, mother)

> family

[1] a.e a.e a.e a.e b.f b.f b.f b.f c.g c.g c.g c.g d.h d.h d.h d.h

Levels: a.e b.e c.e d.e a.f b.f c.f d.f a.g b.g c.g d.g a.h b.h c.h d.h

> family <- interaction(father, mother)

> fit <- lm.rrpp(y~ A / family + B / family + A * B * C, data = df)

Warning: Because variables in the linear model are redundant,

the linear model design has been truncated (via QR decomposition).

Original X columns: 17

Final X columns (rank): 8

Check coefficients or degrees of freedom in ANOVA to see changes.

> anova(fit)

Warning ANOVA is missing some terms, likely because

some independent variables were redundant.

If the residual SS is 0, results should not be trusted

Analysis of Variance, using Residual Randomization

Permutation procedure: Randomization of null model residuals

Number of permutations: 1000

Estimation method: Ordinary Least Squares

Sums of Squares and Cross-products: Type I

Effect sizes (Z) based on F distributions

Df SS MS Rsq F Z Pr(>F)

A 1 0.3134 0.31339 0.05061 1.1107 0.53677 0.298

B 1 0.9008 0.90082 0.14549 3.1926 1.22115 0.114

C 1 0.0057 0.00567 0.00091 0.0201 -1.25962 0.886

A:family 1 1.6208 1.62083 0.26178 5.7443 1.71771 0.034 *

A:C 1 1.0806 1.08057 0.17452 3.8296 1.28432 0.095 .

B:C 1 0.0017 0.00168 0.00027 0.0060 -1.66583 0.948

A:B:C 1 0.0113 0.01133 0.00183 0.0402 -1.04955 0.843

Residuals 8 2.2573 0.28216 0.36458

Total 15 6.1916

---

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Call: lm.rrpp(f1 = y ~ A/family + B/family + A * B * C, data = df)

Sure enough, the same ANOVA but the terms have changed. A:B has disappeared because there is no A:B interaction without replication. One might wish to add family to the model. Not a problem.

> fit <- lm.rrpp(y~ family + A / family + B / family + A * B * C, data = df)

Warning: Because variables in the linear model are redundant,

the linear model design has been truncated (via QR decomposition).

Original X columns: 17

Final X columns (rank): 8

Check coefficients or degrees of freedom in ANOVA to see changes.

> anova(fit)

Warning ANOVA is missing some terms, likely because

some independent variables were redundant.

If the residual SS is 0, results should not be trusted

Analysis of Variance, using Residual Randomization

Permutation procedure: Randomization of null model residuals

Number of permutations: 1000

Estimation method: Ordinary Least Squares

Sums of Squares and Cross-products: Type I

Effect sizes (Z) based on F distributions

Df SS MS Rsq F Z Pr(>F)

family 3 2.8350 0.94501 0.45789 3.3492 1.4530 0.072 .

C 1 0.0057 0.00567 0.00091 0.0201 -1.2191 0.877

C:A 1 1.0806 1.08057 0.17452 3.8296 1.2843 0.095 .

C:B 1 0.0017 0.00168 0.00027 0.0060 -1.6658 0.948

C:A:B 1 0.0113 0.01133 0.00183 0.0402 -1.0495 0.843

Residuals 8 2.2573 0.28216 0.36458

Total 15 6.1916

---

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Call: lm.rrpp(f1 = y ~ family + A/family + B/family + A * B * C, data = df)

>

With more data, maybe the ANOVA would look different, but This illustrates that if family is a variable in the model rather than only nested within A or B, A and B become redundant. This might seem awkward, but notice we get the same answer with:

> fit <- lm.rrpp(y~ family + (A * B * C) / family, data = df)

Warning: Because variables in the linear model are redundant,

the linear model design has been truncated (via QR decomposition).

Original X columns: 35

Final X columns (rank): 8

Check coefficients or degrees of freedom in ANOVA to see changes.

> anova(fit)

Warning ANOVA is missing some terms, likely because

some independent variables were redundant.

If the residual SS is 0, results should not be trusted

Analysis of Variance, using Residual Randomization

Permutation procedure: Randomization of null model residuals

Number of permutations: 1000

Estimation method: Ordinary Least Squares

Sums of Squares and Cross-products: Type I

Effect sizes (Z) based on F distributions

Df SS MS Rsq F Z Pr(>F)

family 3 2.8350 0.94501 0.45789 3.3492 1.4530 0.072 .

C 1 0.0057 0.00567 0.00091 0.0201 -1.2191 0.877

C:A 1 1.0806 1.08057 0.17452 3.8296 1.2843 0.095 .

C:B 1 0.0017 0.00168 0.00027 0.0060 -1.6658 0.948

C:A:B 1 0.0113 0.01133 0.00183 0.0402 -1.0495 0.843

Residuals 8 2.2573 0.28216 0.36458

Total 15 6.1916

---

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Call: lm.rrpp(f1 = y ~ family + (A * B * C)/family, data = df)

If it were me, this is where I might settle. A is the environment of one parent. B is the environment of the other parent. Males and females come from both environments but pairs cannot fully cross parental environments. C is the environment where offspring are raised. We have the main effect C, which is what we are interested in, experimentally. We have all factor interactions, and we know that family is nested within the interactions.

This is not my area of study so there might be a better argument for a different model. Regardless, when nesting random effects, it seems the design should work, but yes, parameters will be lost from the model for interactions that cannot exist.

Hope this is helpful.

Mike

To view this discussion on the web, visit https://groups.google.com/d/msgid/geomorph-r-package/735d0458-8a94-44be-8992-5939806fbec0n%40googlegroups.com.

Tyler Serdar

unread,

Mar 3, 2022, 8:15:37 AM3/3/22

to geomorph R package

Dear Mike,

thank you very much for your extensive explanation and the many examples, which improve my understanding of how to model random effects in this function. However, this is not a quantitative genetics study, and I feel that I have failed to explain the experiment properly, seeing your data simulation, model design and your statement that we are only interested in the main effect C. However, the interest of the study is whether shape variation in individual offspring is best explained by maternal, paternal or offspring environment (or their interactions) - while controlling for family identity. We have actively manipulated all maternal, paternal and offspring environments.

I have simulated and attached an example dataset to this post. The model here would simply be:

y ~ maternal.environment * paternal environment * offspring.environment.

Now, the challenge lies in properly accounting for family effects so as to understand the independent and combined effects of these three fixed effects. Unfortunately, your example focuses only on the offspring environment - but we are interested in whether variation in maternal and paternal environments likewise explain offspring shape.

A similar analysis, using a normal mixed-effect linear models is presented in Table 1 of this paper: https://doi.org/10.3389/fevo.2020.591074

Sorry for my failure to explain the experimental set-up aedaquately before.

Cheers,

Tyler

example.csv

Mike Collyer

unread,

Mar 3, 2022, 9:41:58 AM3/3/22

to geomorph R package

OK, with your data, this seemed to work for me:

> rdf <- rrpp.data.frame(A= df$mother.environment, B = df$father.environment, C = df$offspring.environment, family = df$family, y = df$y)

> fit <- lm.rrpp(y~(A * B * C) / family + family, data = rdf, SS.type = "III")

Warning: Because variables in the linear model are redundant,

the linear model design has been truncated (via QR decomposition).

Original X columns: 98

Final X columns (rank): 22

Check coefficients or degrees of freedom in ANOVA to see changes.

> anova(fit)

Warning ANOVA is missing some terms, likely because

some independent variables were redundant.

If the residual SS is 0, results should not be trusted

Analysis of Variance, using Residual Randomization

Permutation procedure: Randomization of null model residuals

Number of permutations: 1000

Estimation method: Ordinary Least Squares

Sums of Squares and Cross-products: Type III

Effect sizes (Z) based on F distributions

Df SS MS Rsq F Z Pr(>F)

A 1 0.1351 0.135105 0.01629 1.4829 0.80509 0.229

B 1 0.0141 0.014148 0.00171 0.1553 -0.48959 0.682

C 1 0.0908 0.090832 0.01095 0.9970 0.49621 0.332

family 8 0.3261 0.040767 0.03932 0.4475 -1.28925 0.904

A:C 1 0.0336 0.033568 0.00405 0.3684 -0.13893 0.550

B:C 1 0.0005 0.000525 0.00006 0.0058 -1.63119 0.936

A:B:C 1 0.0001 0.000107 0.00001 0.0012 -1.99598 0.976

A:B:C:family 7 0.2940 0.042003 0.03545 0.4610 -1.03453 0.854

Residuals 71 6.4685 0.091106 0.77987

Total 92 8.2943

Call: lm.rrpp(f1 = y ~ (A * B * C)/family + family, SS.type = "III", data = rdf)

Maybe this is what you are looking to do? Note, you could try to force more interactions, but it does not work:

> fit <- lm.rrpp(y~(A * B) / family + (A* C) / family + (B * C) / family + family, data = rdf, SS.type = "III")

Warning: Because variables in the linear model are redundant,

the linear model design has been truncated (via QR decomposition).

Original X columns: 87

Final X columns (rank): 22

Check coefficients or degrees of freedom in ANOVA to see changes.

> anova(fit)

Warning ANOVA is missing some terms, likely because

some independent variables were redundant.

If the residual SS is 0, results should not be trusted

Analysis of Variance, using Residual Randomization

Permutation procedure: Randomization of null model residuals

Number of permutations: 1000

Estimation method: Ordinary Least Squares

Sums of Squares and Cross-products: Type III

Effect sizes (Z) based on F distributions

Df SS MS Rsq F Z Pr(>F)

A 1 0.1351 0.135105 0.01629 1.4829 0.80509 0.229

B 1 0.0141 0.014148 0.00171 0.1553 -0.48959 0.682

C 1 0.0908 0.090832 0.01095 0.9970 0.49621 0.332

family 8 0.4092 0.051156 0.04934 0.5615 -0.85688 0.802

A:C 1 0.0336 0.033568 0.00405 0.3684 -0.13893 0.550

B:C 1 0.0005 0.000525 0.00006 0.0058 -1.63119 0.936

A:C:family 8 0.2941 0.036762 0.03546 0.4035 -1.36748 0.905

Residuals 71 6.4685 0.091106 0.77987

Total 92 8.2943

Call: lm.rrpp(f1 = y ~ (A * B)/family + (A * C)/family + (B * C)/family + family, SS.type = "III", data = rdf)

>

The reason for this is the many columsn of 0s in the model design matrix, i.e.,

> df2 <- as.data.frame(rdf[c("A", "B", "family")])

> model.matrix(~(A * B) / family, data = df2)

(Intercept) AB BB AB:BB AA:BA:familyf10 AB:BA:familyf10 AA:BB:familyf10 AB:BB:familyf10 AA:BA:familyf11

1 1 0 0 0 0 0 0 0 0

2 1 0 0 0 0 0 0 0 0

3 1 0 0 0 0 0 0 0 0

4 1 0 0 0 0 0 0 0 0

5 1 0 0 0 0 0 0 0 0

6 1 0 0 0 0 0 0 0 0

7 1 0 0 0 0 0 0 0 0

8 1 0 0 0 0 0 0 0 0

9 1 0 0 0 0 0 0 0 0

10 1 0 0 0 0 0 0 0 0

11 1 0 0 0 0 0 0 0 0

12 1 0 0 0 0 0 0 0 0

13 1 0 0 0 0 0 0 0 0

14 1 0 0 0 0 0 0 0 0

15 1 0 0 0 0 0 0 0 0

16 1 1 1 1 0 0 0 0 0

17 1 1 1 1 0 0 0 0 0

18 1 1 1 1 0 0 0 0 0

19 1 1 1 1 0 0 0 0 0

20 1 1 1 1 0 0 0 0 0

21 1 1 1 1 0 0 0 0 0

22 1 1 1 1 0 0 0 0 0

AB:BA:familyf11 AA:BB:familyf11 AB:BB:familyf11 AA:BA:familyf2 AB:BA:familyf2 AA:BB:familyf2

1 0 0 0 0 0 0

2 0 0 0 0 0 0

3 0 0 0 0 0 0

4 0 0 0 0 0 0

5 0 0 0 1 0 0

6 0 0 0 1 0 0

7 0 0 0 0 0 0

8 0 0 0 0 0 0

9 0 0 0 1 0 0

10 0 0 0 1 0 0

11 0 0 0 1 0 0

12 0 0 0 1 0 0

13 0 0 0 1 0 0

14 0 0 0 1 0 0

15 0 0 0 1 0 0

16 0 0 0 0 0 0

17 0 0 0 0 0 0

18 0 0 0 0 0 0

19 0 0 0 0 0 0

20 0 0 0 0 0 0

21 0 0 0 0 0 0

22 0 0 0 0 0 0

AB:BB:familyf2 AA:BA:familyf3 AB:BA:familyf3 AA:BB:familyf3 AB:BB:familyf3 AA:BA:familyf4 AB:BA:familyf4

1 0 0 0 0 0 0 0

2 0 0 0 0 0 0 0

3 0 0 0 0 0 0 0

4 0 0 0 0 0 0 0

5 0 0 0 0 0 0 0

6 0 0 0 0 0 0 0

7 0 0 0 0 0 0 0

8 0 0 0 0 0 0 0

9 0 0 0 0 0 0 0

10 0 0 0 0 0 0 0

11 0 0 0 0 0 0 0

12 0 0 0 0 0 0 0

13 0 0 0 0 0 0 0

14 0 0 0 0 0 0 0

15 0 0 0 0 0 0 0

16 0 0 0 0 1 0 0

17 0 0 0 0 1 0 0

18 0 0 0 0 1 0 0

19 0 0 0 0 0 0 0

20 0 0 0 0 0 0 0

21 0 0 0 0 0 0 0

22 0 0 0 0 1 0 0

AA:BB:familyf4 AB:BB:familyf4 AA:BA:familyf5 AB:BA:familyf5 AA:BB:familyf5 AB:BB:familyf5 AA:BA:familyf6

1 0 0 0 0 0 0 0

2 0 0 0 0 0 0 0

3 0 0 0 0 0 0 0

4 0 0 0 0 0 0 0

5 0 0 0 0 0 0 0

6 0 0 0 0 0 0 0

7 0 0 0 0 0 0 0

8 0 0 0 0 0 0 0

9 0 0 0 0 0 0 0

10 0 0 0 0 0 0 0

11 0 0 0 0 0 0 0

12 0 0 0 0 0 0 0

13 0 0 0 0 0 0 0

14 0 0 0 0 0 0 0

15 0 0 0 0 0 0 0

16 0 0 0 0 0 0 0

17 0 0 0 0 0 0 0

18 0 0 0 0 0 0 0

19 0 1 0 0 0 0 0

20 0 1 0 0 0 0 0

21 0 1 0 0 0 0 0

22 0 0 0 0 0 0 0

AB:BA:familyf6 AA:BB:familyf6 AB:BB:familyf6 AA:BA:familyf7 AB:BA:familyf7 AA:BB:familyf7 AB:BB:familyf7

1 0 0 0 0 0 0 0

2 0 0 0 0 0 0 0

3 0 0 0 0 0 0 0

4 0 0 0 0 0 0 0

5 0 0 0 0 0 0 0

6 0 0 0 0 0 0 0

7 0 0 0 0 0 0 0

8 0 0 0 0 0 0 0

9 0 0 0 0 0 0 0

10 0 0 0 0 0 0 0

11 0 0 0 0 0 0 0

12 0 0 0 0 0 0 0

13 0 0 0 0 0 0 0

14 0 0 0 0 0 0 0

15 0 0 0 0 0 0 0

16 0 0 0 0 0 0 0

17 0 0 0 0 0 0 0

18 0 0 0 0 0 0 0

19 0 0 0 0 0 0 0

20 0 0 0 0 0 0 0

21 0 0 0 0 0 0 0

22 0 0 0 0 0 0 0

AA:BA:familyf8 AB:BA:familyf8 AA:BB:familyf8 AB:BB:familyf8 AA:BA:familyf9 AB:BA:familyf9 AA:BB:familyf9

1 0 0 0 0 0 0 0

2 0 0 0 0 0 0 0

3 0 0 0 0 0 0 0

4 0 0 0 0 0 0 0

5 0 0 0 0 0 0 0

6 0 0 0 0 0 0 0

7 0 0 0 0 0 0 0

8 0 0 0 0 0 0 0

9 0 0 0 0 0 0 0

10 0 0 0 0 0 0 0

11 0 0 0 0 0 0 0

12 0 0 0 0 0 0 0

13 0 0 0 0 0 0 0

14 0 0 0 0 0 0 0

15 0 0 0 0 0 0 0

16 0 0 0 0 0 0 0

17 0 0 0 0 0 0 0

18 0 0 0 0 0 0 0

19 0 0 0 0 0 0 0

20 0 0 0 0 0 0 0

21 0 0 0 0 0 0 0

22 0 0 0 0 0 0 0

AB:BB:familyf9

1 0

2 0

3 0

4 0

5 0

6 0

7 0

8 0

9 0

10 0

11 0

12 0

13 0

14 0

15 0

16 0

17 0

18 0

19 0

20 0

21 0

22 0

[ reached getOption("max.print") -- omitted 71 rows ]

attr(,"assign")

[1] 0 1 2 3 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4

attr(,"contrasts")

attr(,"contrasts")$A

[1] "contr.treatment"

attr(,"contrasts")$B

[1] "contr.treatment"

attr(,"contrasts")$family

[1] "contr.treatment"

>

That’s the best I can suggest, to use:

fit <- lm.rrpp(y~(A * B * C) / family + family, data = rdf, SS.type = "III”)

Best,

Mike

To view this discussion on the web, visit https://groups.google.com/d/msgid/geomorph-r-package/1741d063-256b-4993-b1d4-f1f9b3aac27dn%40googlegroups.com.
<example.csv>

Tyler Serdar

unread,

Mar 3, 2022, 9:47:05 AM3/3/22

to geomorph R package

Thanks for your advice, Mike.

In conclusion, it seems to be impossible to correct for family but also get the A:B interaction (i.e. not to truncate the model) at the same time.

Thus, it appears that I have to tell the reviewer that doing so is simply impossible.