1
$\begingroup$

I compare mean hounsfield units (HU) in two different levels of the spine (v1, v2) at two different time points (time1, time2) in 50 patients in R.

test data:

d1 <- data.frame(patient=c(1,2,3),v1_time1=c(50,50,50), v1_time2=c(48,48,48),v2_time1=c(45,45,45), v2_time2=c(30,29,27))
d1

    patient v1_time1 v1_time2 v2_time1 v2_time2
1       1       50       48       45       30
2       2       50       48       45       29
3       3       50       48       45       27

the null hypothesis would be that there is no difference in v1_time1/v1time2 compared to v2_time1/v2_time2

the alternative hypothesis would be that the ratios differ. (e.g. there is no change in v1, but in v2).

What would be the appropriate test?

I thought about an Chi-squared test or a z-Test (prop.test), however, I am not quite sure if my approach is valid.

Would a simple t.test be appropriate?

t.test(d1$v2_time1,d1$v2_time2)
t.test(d1$v1_time1,d1$v1_time2) 
$\endgroup$
1
  • $\begingroup$ You have a continuous variable (the ratio) to compre according to a cathegorical variable (factor v1 or v2), so, depending on the normal distribution of the ratio (normal distributed or not) I would use the T-test or the Wilcoxon test. $\endgroup$ – R18 Dec 11 '17 at 11:40
1
$\begingroup$

No, I think a t-test won't be appropriate in this case. The main problem is that your HU ratios are not independent measurements, i.e. you measure each patient twice, once for position and once for time. So you'd have to account for this. Another problem may be the fact that you have proportions. However, if your values fall around the center and are not predominantly distributed around the tails of the distribution, i.e. zero and one, you might get away with using the normal distribution to analyse your data.

Here's what I would do given my quick reading of the question:

1.) Restructure your data:

library(tidyverse)
d1 %>%
gather(value="HU", key="Variable", -patient)  -> new_d1

> new_d1
   patient Variable HU
1        1 v1_time1 50
2        2 v1_time1 50
3        3 v1_time1 50
4        1 v1_time2 48
5        2 v1_time2 48
6        3 v1_time2 48
7        1 v2_time1 45
8        2 v2_time1 45
9        3 v2_time1 45
10       1 v2_time2 30
11       2 v2_time2 29
12       3 v2_time2 27

2.) Fit a linear-mixed effects model with lmer() from the lme4 package to account for the dependent measurements in your experiment

library(lme4)
fit <- lmer(HU ~ Variable + (1|patient), new_d1)
# inspect results
summary(fit)

Now at this point, you also want to check whether the variance of the residuals in your model is roughly equal (note: with three points, i.e. patients, this won't show you much right now).

plot(fit)

3.)To test for difference in your estimates, you can use the emmeans() function in the emmeans package:

library(emmeans)
emmeans(fit, pairwise ~ Variable)

> emmeans(m, pairwise~key)
$emmeans
 key        emmean        SE df lower.CL upper.CL
 v1_time1 50.00000 0.4409586  8 48.98315 51.01685
 v1_time2 48.00000 0.4409586  8 46.98315 49.01685
 v2_time1 45.00000 0.4409586  8 43.98315 46.01685
 v2_time2 28.66667 0.4409586  8 27.64981 29.68352

Degrees-of-freedom method: kenward-roger 
Confidence level used: 0.95 

$contrasts
 contrast            estimate        SE df t.ratio p.value
 v1_time1 - v1_time2  2.00000 0.6236095  6   3.207  0.0673
 v1_time1 - v2_time1  5.00000 0.6236095  6   8.018  0.0008
 v1_time1 - v2_time2 21.33333 0.6236095  6  34.209  <.0001
 v1_time2 - v2_time1  3.00000 0.6236095  6   4.811  0.0117
 v1_time2 - v2_time2 19.33333 0.6236095  6  31.002  <.0001
 v2_time1 - v2_time2 16.33333 0.6236095  6  26.192  <.0001

P value adjustment: tukey method for comparing a family of 4 estimates

NOTE: There might be a better or more appropriate way to do that especially when seeing the results using the full data set, but I think the analysis should head in this direction given the dependent structure in your data.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.