На основании предложения пользователя я разместил это на CrossValidated: https://stats.stackexchange.com/questions/468693/r-regression-analysis-between-two-numeric-variables-stratified-by-groups
У меня есть следующий фрейм данных из 48 образцов со следующими столбцами: sample (образец идентификатор), gene (log2 экспрессия гена случайно выбранного гена), label (тренированные и не тренированные мыши), strain (линия мышей), weight (масса мышей), running.time (сколько времени потрачено на беговую дорожку) и VO2max (максимальное потребление O2).
> head(test)
sample label strain weight VO2max running.time gene
1 A81 exercised ANT1 ME 28.8 6192.907 21.00000 18.26548
2 A13 exercised ANT1 ME 26.8 6598.627 22.25000 17.55368
3 A43 exercised ANT1 ME 26.4 6368.918 20.33333 17.45001
4 A76 exercised ANT1 ME 23.0 6947.636 18.41667 17.87199
5 A17 exercised ANT1 ME 25.8 6526.043 21.58333 17.64396
6 A30 exercised ANT1 ME 29.2 6562.106 20.58333 17.81958
> dput(test)
structure(list(sample = c("A81", "A13", "A43", "A76", "A17",
"A30", "A58", "A75", "A82", "A88", "A4", "A18", "I50", "I94",
"I75", "I59", "I63", "I84", "I13", "I95", "I18", "I100", "I62",
"I85", "B69", "B80", "B11", "B9", "B98", "B4", "B70", "B19",
"B85", "B3", "B10", "B47", "E56", "E76", "E12", "E78", "E1",
"E50", "E42", "E64", "E98", "E100", "E7", "E53"), label = c("exercised",
"exercised", "exercised", "exercised", "exercised", "exercised",
"non_exercised", "non_exercised", "non_exercised", "non_exercised",
"non_exercised", "non_exercised", "exercised", "exercised", "exercised",
"exercised", "exercised", "exercised", "non_exercised", "non_exercised",
"non_exercised", "non_exercised", "non_exercised", "non_exercised",
"exercised", "exercised", "exercised", "exercised", "exercised",
"exercised", "non_exercised", "non_exercised", "non_exercised",
"non_exercised", "non_exercised", "non_exercised", "exercised",
"exercised", "exercised", "exercised", "exercised", "exercised",
"non_exercised", "non_exercised", "non_exercised", "non_exercised",
"non_exercised", "non_exercised"), strain = c("ANT1 ME", "ANT1 ME",
"ANT1 ME", "ANT1 ME", "ANT1 ME", "ANT1 ME", "ANT1 ME", "ANT1 ME",
"ANT1 ME", "ANT1 ME", "ANT1 ME", "ANT1 ME", "IAI ME", "IAI ME",
"IAI ME", "IAI ME", "IAI ME", "IAI ME", "IAI ME", "IAI ME", "IAI ME",
"IAI ME", "IAI ME", "IAI ME", "B6 ME", "B6 ME", "B6 ME", "B6 ME",
"B6 ME", "B6 ME", "B6 ME", "B6 ME", "B6 ME", "B6 ME", "B6 ME",
"B6 ME", "EC77 ME", "EC77 ME", "EC77 ME", "EC77 ME", "EC77 ME",
"EC77 ME", "EC77 ME", "EC77 ME", "EC77 ME", "EC77 ME", "EC77 ME",
"EC77 ME"), weight = c(28.8, 26.8, 26.4, 23, 25.8, 29.2, 27.6,
25, 28.6, 29.4, 29, 29.2, 31.6, 29.9, 34, 32.5, 31.7, 31, 28.8,
29.8, 27.9, 32.8, 32.1, 29, 32.2, 28.4, 30, 26.3, 27.8, 30.8,
32.1, 28, 30, 35.3, 29.9, 29.6, 30, 31, 24, 25.8, 27.4, 30, 31,
30.4, 28, 27.7, 26.25, 29.1), VO2max = c(6192.907, 6598.627,
6368.918, 6947.636, 6526.043, 6562.106, 6626.237, 6464.229, 6692.121,
6990.765, 6586.469, 6773.607, 7755.965, 8504.163, 8081.548, 7500.55,
7693.678, 7909.825, 7588.132, 6710.412, 6846.685, 7498.796, 7196.46,
7856.87, 7605.042, 8616.026, 7942.862, 8306.043, 7949.328, 8376.683,
8480.413, 7422.928, 8006.701, 7987.724, 8091.422, 7622.233, 7396.47,
7538.512, 7144.345, 7032.685, 6621.896, 7800.92, 6900.873, 7078.454,
6846.685, 6165.303, 7327.488, 6787.016), running.time = c(21,
22.25, 20.3333333333333, 18.4166666666667, 21.5833333333333,
20.5833333333333, 18.75, 17.5833333333333, 20.9166666666667,
20.9166666666667, 19.3333333333333, 19.6666666666667, 36.75,
35.5, 33.5, 33.9166666666667, 31.3333333333333, 36, 33.0833333333333,
33.1666666666667, 32.4166666666667, 34.1666666666667, 32.1666666666667,
32.1666666666667, 37.3333333333333, 36.1666666666667, 33.6666666666667,
36.4166666666667, 35.4166666666667, 35.5, 36.5833333333333, 33.8333333333333,
34.4166666666667, 33.5833333333333, 35.6666666666667, 34.75,
36.1666666666667, 31.1666666666667, 37.5, 31.6666666666667, 33.1666666666667,
32.0833333333333, 33.5833333333333, 32.1666666666667, 32.4166666666667,
31.25, 33.5833333333333, 31.5), gene = c(18.265481735128, 17.553682163378,
17.450005832624, 17.8719936539426, 17.6439644306893, 17.8195797128757,
17.7327533149751, 17.8626472422709, 18.4375411245124, 17.7206616316039,
17.8324727651774, 17.6831637379254, 16.3797370363944, 16.2274547911749,
16.1081681682918, 15.7581300013604, 16.3283156117628, 16.6112399890101,
16.2097766576023, 16.353875382536, 16.1630593358039, 16.1941586675123,
16.270535309697, 16.4154082664562, 16.0543656477653, 16.2056037427263,
16.0414009536146, 15.8233199504989, 15.9644539321719, 16.0035065577481,
15.1967778769911, 15.5873577174736, 15.4748920755074, 15.8414050703223,
16.036749175609, 15.6950267375145, 16.2487211678851, 16.4280535060295,
16.1600327464311, 15.9933079171589, 16.4159279247705, 15.8740001632907,
16.2547699681823, 15.9199201254046, 15.7748035384056, 16.0895632695324,
15.9950306563657, 16.1150404336735)), class = "data.frame", row.names = c(NA,
-48L))
Я хочу посмотреть, коррелирует ли экспрессия гена gene (numeri c) с непрерывным переменная VO2max (numeri c) для двух групп мышей с упражнениями против non_exercised (label) и четырех линий мышей (strain).
Это модель, о которой я думаю, но не конечно, если я прав:
# fit model
> fit <- lm(gene ~ VO2max + label + strain, data = test)
# output
> tidy(fit)
# A tibble: 6 x 5
term estimate std.error statistic p.value
<chr> <dbl> <dbl> <dbl> <dbl>
1 (Intercept) 18.2 0.686 26.5 7.67e-28
2 VO2max -0.0000434 0.000101 -0.430 6.69e- 1
3 labelnon_exercised -0.110 0.0767 -1.44 1.57e- 1
4 strainB6 ME -1.93 0.176 -11.0 6.16e-14
5 strainEC77 ME -1.70 0.111 -15.3 8.73e-19
6 strainIAI ME -1.53 0.142 -10.8 1.21e-13
Я не уверен, правильный ли это подход или как интерпретировать результаты. Может кто-нибудь помочь?
Спасибо