LSM calculation in 2x2 2-stage design

Phoenix Certara Forums Phoenix WNL basics
1

Hello all!

Could someone kindly explain how LSM is computing in PWNL?

I am use this model:

Sequence+Stage+Period(Stage)+Formulation+SequenceStage+Subject(SequenceStage)

All effects is “fixed”.

And it is clear for me how coefficients of least squares means is computed, but not clear how to estimate LSM.

If I try to calculate LSM like = (Test_mean_seqence1_stage1 + Test_mean_seqence2_stage1 + Test_mean_seqence1_stage2 + Test_mean_seqence1_stage2)/4 - a have result not equal to PWNL.

Is there any analytical expression to do it?

2

Hi Vladimir,

Sorry, there’s no simple formula as for conventional model due to interaction and nested terms

Citing R help:

A special type of linear estimates is the so called least–squares means (or LS–means). Other names for these estimates include population means and marginal means. Consider an imaginary field experiment analyzed with model of the form

lm( y ~ treat + block + year)

where treat is a treatment factor, block is a blocking factor and year is the year (a factor) where the experiment is repeated over several years. This model specifies the conditional mean E(Y |treat, block, year). One may be interested in predictions of the form E(Y |treat). This quantity can not formally be found from the model. However, it is tempting to average the fitted values of E(Y |treat, block, year) across the levels of block and year and think of this average as E(Y |treat). This average is precisely what is called the LS–means. If the experiment is balanced then this average is identical to the average of the observations when stratified according to treat.

So LSmeans are conditional means over other factors. And it’s true, it is just an average in case of simple model.

The problem is that 2-stage EMA model is not simple and you need to take into account intrablock information. No analytical ways out

hope it helps,

Mittyright

3

Hi again

If you are interested in cross-validation with R and in understanding of algo, please review the code below

library(dplyr)
library(nlme)
Subj=as.factor(rep(1:20, each = 2))
Seq=as.factor(c(rep(c("TR","RT"), each = 6), 
                rep(c("RT","TR"), each = 6), 
                rep(c("TR","RT"), each = 4), 
                rep(c("RT","TR"), each = 4)))
Per=as.factor((rep(c(1,2), 20)))
Trt=as.factor(c(rep(c(1,2), 6), 
                rep(c(2,1), 6), 
                rep(c(1,2), 4), 
                rep(c(2,1), 4)))
Stg=as.factor(c(rep(1,24), 
                rep(2,16)))
PK=c(4.8, 4.6, 1.4, 1.3, 3.9, 4.4, 1.5, 2.4, 3.8, 2.9, 3.3, 2.5, 1.7, 1.5, 1.3, 1.9, 2.1, 2.4, 11.4, 13.9, 1.5, 2.3, 2.6, 2.4,
     4.0,4.2,0.7, 0.5, 10.3, 11.9,4.7, 4.8, 5.9,5.5, 3.8, 5.4, 4.3, 2.1, 3.6, 2.4)
M=lm(PK ~ Seq + Stg/Per + Stg:Seq/Subj + Trt)
newdat <- expand.grid(Trt = levels(Trt), Seq = levels(Seq), Subj = levels(Subj), Per = levels(Per), Stg = levels(Stg))
preddata <- cbind(newdat, predict(M, newdat))
dta <- cbind.data.frame(Trt, Seq, Subj, Per,  Stg, PK)
preddatawoNA <- na.omit(left_join(preddata, dta))
preddatawoNA%>%
  group_by(Trt, Seq, Per, Stg) %>%
  summarize(Subjmean = mean(`predict(M, newdat)`)) %>%
  group_by(Trt) %>%
  summarize(lsmean = mean(Subjmean))

All the best,

Mittyright