Я уже искал на Форуме часы (действительно) и начинаю чувствовать слабое чувство, что постепенно схожу с ума, тем более, что мне кажется, что это действительно легко решаемая проблема.
Что делатьЯ хочу сделать?
По сути, я хочу смоделировать клинические данные.В частности, для каждого пациента (столбец 1: ID) произвольная оценка (столбец 3: оценка), зависящая от назначенной группы лечения (столбец 2: группа).
set.seed(123)
# Number of subjects in study
n_patients = 1000
# Score: Mean and SDs
mean_verum = 70
sd_verum = 20
mean_placebo = 40
sd_placebo = 20
# Allocating to Treatment groups:
data = data.frame(id = as.character(1:n_patients))
data$group[1:(n_patients/2)] <- "placebo"
data$group[(n_patients/2+1):n_patients] <- "verum"
# Attach Score for each treatment group
data$score <- ifelse(data$group == "verum", rnorm(n=100, mean=mean_verum, sd=sd_verum), rnorm(n=100, mean=mean_placebo, sd=sd_placebo))
Пока все просто.Теперь я хочу 1) рассчитать вероятность наступления события (функция логита) в зависимости от оценки.Затем, 2) Я хочу фактически назначить Событие, в зависимости от вероятности (rbinom).
Я хочу сделать это для n различных вероятностей / Событий.Это код, который я использовал до сих пор:
Рассчитать вероятности:
a = -1
b = 0.01
p1 = 1-exp(a+b*data$score)/(1+exp(a+b*data$score))
data$p_AE1 <- p1
a = -0.5
b = 0.01
p1 = 1-exp(a+b*data$score)/(1+exp(a+b*data$score))
data$p_AE2 <- p1
…
Назначить события:
data$Abbruch_AE1 <- rbinom(n_patients, 1, data$p_E1)
data$Abbruch_AE2 <- rbinom(n_patients, 1, data$p_E2)
…
Очевидно, что это действительно неэффективно, так какЯ хотел бы легко увеличить или уменьшить это, в зависимости от того, сколько вероятностей / событий я хочу смоделировать.
Проблема в том, что я просто не понимаю, как я могу одновременно а) генерировать новый, единственный столбецв кадре данных, где я хочу ввести значения для каждого, b) выполнить функцию для назначения вероятностей / событий и c) сделать это для числа n различных формул, которые имеют свои конкретные a и b.
Я уверен, что решение этой проблемы простое - мне не удалось сделать все эти вещи одновременно, и я бы хотел, чтобы это было в конечном итоге.Я поиграл с циклами, но все безрезультатно.
Любая помощь будет принята с благодарностью!
Вот так выглядит мой фрейм данных:
structure(list(id = structure(1:3, .Label = c("1", "2", "3"), class = "factor"),
group = c("placebo", "placebo", "placebo"), score = c(25.791868726014,
45.1376741831306, 35.0661624307525), p_AE1 = c(0.677450814266315,
0.633816117436442, 0.656861351663365), p_AE2 = c(0.560226492151216,
0.512153420188678, 0.537265362130761), p_AE3 = c(0.435875409622676,
0.389033483248856, 0.413221988111604), p_AE4 = c(0.319098312196655,
0.278608032377073, 0.299294085148527), p_AE5 = c(0.221332386680766,
0.189789774534235, 0.205762225373345), p_AE6 = c(0.147051201194953,
0.124403316086538, 0.135795233451071), p_AE7 = c(0.0946686004658072,
0.0793379289917946, 0.0870131973838217), p_AE8 = c(0.0596409872667201,
0.0496714832182721, 0.0546471270895262), AbbruchAE1 = c(1L,
1L, 1L), AbbruchAE2 = c(1L, 1L, 0L), AbbruchAE3 = c(0L, 0L,
0L), AbbruchAE4 = c(0L, 1L, 0L), AbbruchAE5 = c(1L, 0L, 0L
), AbbruchAE6 = c(1L, 0L, 0L), AbbruchAE7 = c(0L, 0L, 0L),
AbbruchAE8 = c(0L, 0L, 0L)), .Names = c("id", "group", "score", "p_AE1", "p_AE2", "p_AE3", "p_AE4", "p_AE5", "p_AE6", "p_AE7", "p_AE8", "AbbruchAE1", "AbbruchAE2", "AbbruchAE3", "AbbruchAE4", "AbbruchAE5", "AbbruchAE6", "AbbruchAE7", "AbbruchAE8"), row.names = c(NA, 3L), class = "data.frame")