※ 引述《yummy7922 (crucify)》之銘言：
: ※ 引述《yummy7922 (crucify)》之銘言：
: : [問題類型]:
: : 程式諮詢(我想用R 做某件事情，但是我不知道要怎麼用R 寫出來)
: : [軟體熟悉度]:
: : 入門(寫過其他程式，只是對語法不熟悉)
: : [問題敘述]:
: : 我的資料是一個重複測量的資料，資料中有當期是否使用藥物的資料(若有使用，設
為1
: : 沒有為0)、還有記錄該筆資料為該位病人的第幾筆觀察值，我希望能算出，藥物轉換
的
: : 率，例如有多少人第一期使用a藥物，但在第二期時轉換成使用b藥物。
: : 資料中共有6種藥物，一位病人最多有48筆觀察值，
: : (6x6)種轉換可能 x 47(個時間隔) = 1692 個機率值。
: : 想請教各位高手們，該怎麼做比較有效率。
: 不好意思，我沒有說明清楚，
: 我想算的機率其實是很直觀的條件機率，例如：
: P(第二期使用b藥物|第一期使用a藥物)
: = (第一期使用a藥且第二期使用b藥的人數)/(第一期使用a藥的人數)
: 但是我不知道該怎麼在分組之後，還能夠給定條件，
: 計算出第一期使用a藥第二期使用b藥的人數。
library(data.table)
library(dplyr)
library(magrittr)
library(Matrix)
# data generation
N_patient = 1000
dat = sample(48, N_patient, replace = TRUE) %>% {
cbind(rep(1:N_patient, times=.), sapply(., seq, from = 1) %>% unlist())
} %>% cbind(sample(6, nrow(.), replace = TRUE)) %>%
data.table() %>% setnames(c("id", "obs_time", "dose"))
dat_transform = dat$dose %>%
spMatrix(length(.), 6, 1:length(.), ., rep(1,length(.))) %>%
as.matrix() %>% data.table() %>% setnames(paste0("dose_", 1:6)) %>%
cbind(select(dat, 1:2), .)
# dat_transform為你的資料格式，透過下面可以轉成dat
dat = dat_transform %>%
mutate(dose = dose_1 + dose_2 * 2 + dose_3 * 3 + dose_4 * 4 +
dose_5 * 5 + dose_6 * 6) %>% select(c(1,2, 9))
# 你的資料column順序跟我的不一樣，select要改不同數字
st = proc.time()
dat_count = dat %>% group_by(id) %>%
mutate(previous_dose = c(0, dose[-length(dose)])) %>%
filter(obs_time > 1) %>%
group_by(obs_time, dose, previous_dose) %>%
summarise(count = n())
transitMatrix = dat_count %>% group_by(obs_time, previous_dose) %>%
mutate(transitP = count / sum(count)) %>% ungroup() %>%
split(.$obs_time) %>%
lapply(function(x) spMatrix(6, 6, x$previous_dose, x$dose,
x$transitP))
proc.time() - st
# user system elapsed
# 0.72 0.05 0.78
隨機過程裡面應該有定義這個矩陣：transition probability matrix
翻譯應該是轉移機率矩陣，下面簡稱轉移矩陣
> transitMatrix[[1]] # 第二期的轉移矩陣
[1,] 0.1595092 0.1288344 0.1411043 0.2085890 0.1901840 0.1717791
[2,] 0.1506024 0.2108434 0.1506024 0.1626506 0.2108434 0.1144578
[3,] 0.1676647 0.1137725 0.1616766 0.1497006 0.2035928 0.2035928
[4,] 0.2000000 0.1161290 0.1806452 0.1612903 0.1870968 0.1548387
[5,] 0.1707317 0.2256098 0.1890244 0.1341463 0.1463415 0.1341463
[6,] 0.1296296 0.1851852 0.2160494 0.1543210 0.1728395 0.1419753
(1,1)位置就是上一期服dose 1，下一期也是dose 1的條件機率
(1,2)位置就是上一期服dose 1，下一期服dose 2的條件機率
剩下依此類推
PS: 因此，轉移矩陣的每一個列，其和為1
備註：
此處用spMatrix不是因為output是sparse matrix (如果人少，則是)
而是程式比較好寫，直接透過矩陣位置做initial，省功夫去排列而已
特此說明。