fitted values in mmrm about mmrm HOT 8 CLOSED

Hi @davood1983 , can you please describe the question/problem?

from mmrm.

Hi,
Thanks for getting back to me.
I have a clinical experiment with 6 assessments.
The first assessment is before applying the treatment so the means in the first assessment are equal.
The model that I fit is
mm2=mmrm(l.colony~ f.day+f.day:Treatment+ar1h(f.day|f.tu), data=d4)
and then change the correlation matrix.
This model imposes the equality of the means at assessment 1, however, when I check the fitted(mm2) this condition is not imposed.
The other problem is that when I change the correlation matrix, the fitted values remain constant.
Is mmrm compatible with the fitted(.) command?
If not, what is your suggestion for getting the fitted values?

Best,
Davood

from mmrm.

Hi Davood,

cool thanks for the initial information. Can you please provide a reproducible example, e.g. with fake data and the full code so I can run this and have a look?

The mmrm package supports the fitted method, see https://openpharma.github.io/mmrm/latest-tag/reference/mmrm_tmb_methods.html .
Please also have a look at emmeans to obtain least square means, I am not sure if that is what you are looking for: https://openpharma.github.io/mmrm/latest-tag/articles/introduction.html?q=emmeans#support-for-emmeans

from mmrm.

library(mmrm)

# Set seed for reproducibility
set.seed(123)

# Number of levels
num_levels <- 4  # 3 treatment levels + 1 control

# Number of replicates for each level
num_replicates <- 5

# Number of assessments
num_assessments <- 6

# Mean values for each level at the first assessment
initial_means <- rep(10, num_levels)

# Simulate data
data <- data.frame(
  Subject = rep(1:(num_replicates * num_levels), each = num_assessments),
  Assessment = rep(1:num_assessments, times = num_replicates * num_levels),
  Treatment = rep(1:num_levels, each = num_replicates),
  Value = numeric(num_replicates * num_levels * num_assessments)
)

# Set initial values for each level
data$Value[data$Assessment == 1] <- rep(initial_means, each = num_replicates)

# Introduce variability to simulate treatment effect
treatment_effect <- c(0, 2, -1, 1)  # Adjust these values as needed
for (i in 2:num_assessments) {
  for (j in 1:num_levels) {
    data$Value[data$Assessment == i & data$Treatment == j] <- 
      data$Value[data$Assessment == (i - 1) & data$Treatment == j] + 
      rnorm(num_replicates, mean = treatment_effect[j], sd = 2)
  }
}

# Print the first few rows of the simulated dataset
head(data)
data$Assessment<-as.factor(data$Assessment)
data$Treatment<-as.factor(data$Treatment)
data$Subject<-as.factor(data$Subject)
n <- length((levels(data$Treatment)))

library(randomcoloR)
palette <- distinctColorPalette(n)

interaction.plot(
  x.factor = data$Assessment,
  trace.factor = data$Treatment,
  response = data$Value,
  col = palette,
  type = "b"
)
mm33 <- mmrm(Value~  Assessment+Assessment:Treatment+cs(Assessment|Subject), data=data)
summary(mm33)
data$c1 <- fitted(mm33)

d5 <- data[data$Assessment=="1",]

library(dplyr)
d5 |> 
  group_by(Treatment) |> 
  summarise(mean(c1))


interaction.plot(
  x.factor = data$Assessment,
  trace.factor = data$Treatment,
  response = c1,
  col = palette,
  type = "b"
)

I have a question, when the first level of the assessment is considered as the reference group, why there is an interaction between the reference group of the assessment and the treatment?

I think the reference group gets the zero value and the interactions should not be estimated in output.
in my real data, these interactions are big numbers.

from mmrm.

Thanks @davood1983 , I just ran your example.

One comment: In the doubly nested loop I get 40 warnings:

There were 40 warnings (use warnings() to see them)
> warnings()
Warning messages:
1: In data$Value[data$Assessment == (i - 1) & data$Treatment ==  ... :
  longer object length is not a multiple of shorter object length

so it might well not do you want.

Also I needed to add some packages that were not included, I edited the code accordingly above.
However c1 is not defined so the last statement does not pass.

If you would like a good answer from us, please ensure the example is really reproducible (run it in a fresh R session - does it work?) and please clearly describe the problem that you see or the question that you have. At the moment I still don't see what is the issue.

from mmrm.

Hi @davood1983 - if this is still an issue, please let me know, if I can run your example I am sure we can fix it!

from mmrm.

from mmrm.

Sounds good, thanks @davood1983 !

from mmrm.