<- function(size, grp = 2, T_2_C = "1:1"){
simple_random set.seed(20210412)
<- seq(1, size, 1)
id_num <- runif(n = size,
random_seq min = 0, max = 1)
<- rank(random_seq)
int_rank
<- as.numeric(substr(T_2_C, 1, 1))
ratio_T <- as.numeric(substr(T_2_C, 3, 3))
ratio_C
if (grp == 2) {
<- ifelse(int_rank <=
group /(ratio_T + ratio_C),
size"T", "C")
else if (grp > 3) {
} <- cut(int_rank, breaks = grp,
group labels = paste("Group", 1:grp))
}
<- data.frame("ID" = id_num,
df "RandomNum" = random_seq,
"Rank" = int_rank, "Group" = group)
# write.csv(df, "simple randomization table.csv",
# row.names = FALSE)
return(df)
}
6 番外篇(一):随机化分组
6.1 引子
实验研究作为因果推断效能最高的一种研究设计,被广泛应用于医学及社会科学研究中。在实验研究中又以随机对照试验(Randomised Controlled Trials,RCT)最为常见。RCT研究一个最明显的特点就是随机化(Randomization)施加干预(李立明 2007),随机化和盲法是RCT研究的基石(Kendall 2003)。随机化干预的重要性在于可以有效避免选择或其他偏倚,提供统一的统计分析方法,以及为精确检验显著性和区间估计提供基础(Roberts and Torgerson 1998; Cox 2009)。
在具体的研究中,随机化干预的基础是对干预对象进行随机分组,然而,在现有针对本科生和研究生阶段学习的统计教材中,却缺少对随机分组方法以及其具体实现方法的详细介绍,一定程度上不利于青年研究者将RCT方法应用到具体的研究实践中。尽管已有文献报道了基于SAS和Stata的随机分组方法(胡良平 et al. 2011; 刘玉秀 et al. 2001; 吴春霖, 王镭, and 李卫兵 2013),但是相较于SAS和Stata而言,R语言安装与使用更为简便,并且由于其开源免费的特性更易于获取。因此,本文基于已有文献,详细介绍三种随机分组的具体方法和操作过程,以及提供基于R语言的实现代码,供其他学者在具体研究中作为参考。
6.2 随机分组方法简介
随机化的概念最早由Fisher于1926年提出(Aylmer Fisher 1926),其后被广泛应用于实验研究中,已成为控制干预组和对照组之间混杂偏倚的必要手段。随机化具体指将受试者按照相同的概率分配在干预组和非干预组的过程(Fleiss, Levin, and Paik 2003),根据实验设计的样本大小以及受试者招募时间长短等因素的不同,随机化分组方法很多种,但较为常用的有三种:简单随机分组(simple randomization)、区组随机分组(block randomization)和分层随机分组(stratified randomization)(Kang, Ragan, and Park 2008)。具体如下:
简单随机分组是指按照一组随机序列将受试者以相同的概率分配至不同的实验组中,最典型的简单随机分组方法是抛硬币(flipping a coin)。简单随机分组操作容易,但是存在一个明显的不足,即当受试者样本量较小时,无法保证纳入不同实验组的受试者数量相等(Kang, Ragan, and Park 2008; Kernan et al. 1999)。一般认为当受试者样本量小于200时不宜采用简单随机分组(Kang, Ragan, and Park 2008)。实际上,在具体的研究过程中(如临床药物试验),通常要求干预组和对照组所分配的受试者数量相等或达到一定的比例,因此即使样本量大于200简单随机分组也不适用。然而,在实际研究过程中,会对简单随机分组方法进行一定的变化以满足研究设计需要,具体见下节。
区组随机分组是指先根据受试者样本量大小及实验分组数量和组间比例设定好区组(block)的长度和数量,将特征相似或相同(如入组时间、地区来源)的受试者编入同一区组,然后在每个区组内部进行简单随机分组。区组随机分组可以有效保证不同的实验组中所分配的受试者数量严格满足研究设计要求。不仅如此,区组随机分组的另一个优势在于,当入组后的受试者出现某种不满足实验要求而需退出实验时,可以直接将该受试者所在的区组舍弃另外补充新的区组,而不破坏整个实验的随机化。
分层随机分组是指先按已知的混杂因素对受试者进行分层后,再在各个层内采用简单随机或者区组随机方法将受试者分配至不同的实验组。分层随机分组可以在实验设计阶段就对部分混杂因素进行控制,从而提高不同实验组中受试者的基线均衡性(Kernan et al. 1999)。
其他随机分组方法还有:协变量自适应随机分组(covariate adaptive randomization)和响应自适应随机分组( response-adaptive randomization)。需要补充的是,在实际研究过程中,很难一次性招募足够的受试者,然后再为受试者分配随机号进行随机分组工作。通常情况是会设定一段时间(一个月或三个月,视具体情况而定)的招募入组期,当有受试者满足入组条件后就需要进行随机分组,既而接受相应的干预,因此,这就需要在受试者入组之前准备好随机分组规则,受试者根据进入研究的顺序,依据预先准备好的随机分组表进入相应的实验组。
6.3 随机分组的具体过程及R语言代码
6.3.1 简单随机分组(simple randomization)
简单随机分组的基本步骤一共有五点:1)按照实验设计的受试者样本量生成受试者编号;2)利用随机函数生成与受试者编号数量相等的随机数字;3)对生成的随机数字进行编秩;4)根据实验设计的组数及分组比例,对随机数字的秩序按照大小或者奇偶数进行分组;5)获得符合实验设计要求的随机分组表。本文提供了一段基于R语言的自编函数,能够满足两组及以上分组且不同组间比例的简单随机分组需求,自编函数代码及具体含义如下:
函数名为simple_random,函数共有3个参数,分别为size(样本量大小)、grp(分组数)和T_2_C(干预组与对照组比例),grp和T_2_C为默认参数,分别默认赋值为2和”1:1”;
set.seed()函数用于设定随机种子数,以便下一步利用随机函数生成的随机数可再现,此步骤在利用计算机进行实验研究的随机化分组时十分关键,因为只有固定了随机种子数,才可以保证实验研究中生成的随机分组表可以复现,这一点对于研究论文的发表和临床药物试验质量稽查十分重要。利用计算机生成的随机数本质为伪随机数(Pseudorandomness),它是通过特定的算法计算所得,其满足随机数的统计特征,在相同的计算平台下可以重复出现。
seq()函数用于生成从1开始的受试者编号;
runif()函数用于生成与受试者人数相等的服从均匀分布的0到1之间随机数字,其中min和max参数可以任意指定,此处也可采用R语言中自带的其他随机函数,如rnorm()等,也可以采用sample()函数直接生成非负随机整数;
rank()函数用于对生成的随机数字按照升序进行编秩;
ratio_T和ratio_C变量用于提取T_2_C参数的赋值,用于分组比例的计算;
if条件语句中的内容即是对随机数字的秩进行分组,其规则为:当为两组且组间比例为1:1时,随机数字的秩位于前50%的受试者编号设定为干预组(“T”),后50%对应的受试者编号设定为对照组(“C”),当为两组且组间比例为1:N时,随机数字的秩位于前1/(1+N)的受试者编号设定为干预组(“T”),剩余受试者编号设定为对照组(“C”),当为三组及以上时,本文提供的自编函数仅支持组间比例相等的情况,其分组规则为将随机数字的秩进行等分,与秩对应的受试者编号分别分入”Group 1”、“Group 2”…“Group N”;
data.frame()和write.csv()函数用于生成随机分组表并输出为csv格式,随机分组表中包含受试者编号、随机数字、随机数字的秩、分组。
实例一:为检验某一治疗高血压新药A的临床疗效,选择已上市B药作为对照组,开展单中心双臂临床试验,需要受试者240名,每组各120名受试者,请按简单随机分组方法制定随机分组表。利用本文提供的自编函数simple_random()则可将参数size设定为240, 参数grp设定为2,参数T_2_C设定为”1:1”(注意此处为英文半角引号),具体如下:
<- simple_random(size = 240, grp = 2, T_2_C = "1:1") res
函数输出的随机分组表见 表 tbl-random-1 ,其中分组为T和C的受试者分别为120名。
library(knitr)
kable(res[c(1:12, 229:240),])
ID | RandomNum | Rank | Group | |
---|---|---|---|---|
1 | 1 | 0.8323749 | 198 | C |
2 | 2 | 0.9552218 | 228 | C |
3 | 3 | 0.5978788 | 134 | C |
4 | 4 | 0.3507679 | 75 | T |
5 | 5 | 0.4315742 | 91 | T |
6 | 6 | 0.6332632 | 147 | C |
7 | 7 | 0.7801558 | 189 | C |
8 | 8 | 0.4699095 | 102 | T |
9 | 9 | 0.3853540 | 81 | T |
10 | 10 | 0.6336118 | 148 | C |
11 | 11 | 0.7365508 | 178 | C |
12 | 12 | 0.4967514 | 106 | T |
229 | 229 | 0.6637343 | 156 | C |
230 | 230 | 0.9801347 | 235 | C |
231 | 231 | 0.1659937 | 36 | T |
232 | 232 | 0.4225586 | 90 | T |
233 | 233 | 0.3599470 | 78 | T |
234 | 234 | 0.6065274 | 136 | C |
235 | 235 | 0.5213893 | 111 | T |
236 | 236 | 0.8585871 | 202 | C |
237 | 237 | 0.6808066 | 161 | C |
238 | 238 | 0.2789617 | 63 | T |
239 | 239 | 0.9939375 | 236 | C |
240 | 240 | 0.7866525 | 190 | C |
6.3.2 区组随机分组(block randomization)
6.3.2.1 方法一
区组随机分组一般有两种办法,第一种是先根据实验分组情况设定区组长度和区组类型组合,再根据区组数量对区组类型进行随机抽样,从而获得每个区组内部的分组情况(Kang, Ragan, and Park 2008),其基本步骤为:1)按照实验设计的受试者样本量生成受试者编号;2)设定区组长度和区组组合;3)根据样本量计算区组数;4)根据区组数从区组组合中随机抽取对应数量的区组类型;5)将随机抽取的区组类型顺序连接,从而获得符合实验设计要求的随机分组表。本文同样提供了一段基于R语言的自编函数,但仅能满足组间比例为1:1的区组随机分组需求,自编函数代码及具体含义如下:
<- function(size,
block_random_m1
block_len,
block_category){
<- seq(1, size, 1)
id_num <- size / block_len
block <- c()
block_set <- c()
block_num for (i in 1:block) {
set.seed(20210412+i)
<- sample(block_category, 1,
block replace = FALSE)
<- block
block_set[i] <- append(block_num,
block_num rep(i, block_len))
}<- as.vector(
group unlist(
strsplit(
paste0(block_set,
collapse=";"),
split=";")
)
)<- data.frame("ID" = id_num,
df "Block" = block_num,
"Group" = group)
# write.csv(df,
# "block radomization table 1.csv",
# row.names = FALSE)
return(df)
}
函数名为block_random_m1,函数共有3个参数,分别为size(样本量大小)、block_len(区组长度)和block_category(区组组合);
seq()函数和set.seed()作用同上;
block变量为根据样本量和区组长度计算的区组数量;
两组实验中常用的区组长度为4,那么参数block_ken即为4,对应的区组组合共有6种情况(“TTCC”, “TCTC”, “TCCT”, “CCTT”, “CTCT”, “CTTC”),其中T表示干预组,C表示对照组,用参数block_category来指定;
for循环中即是依据计算所得的区组数量,利用sample()函数从全部区组组合中随机抽取相应数量的区组组合,其中sample()函数的第二个参数用于指定每次抽取的组合数,本文中采取每次抽取1个,第三个参数replace用于指定是否为有放回抽样,当每次抽取数量大于1时应采取无放回抽样,replace的值应设定为FALSE;
block_set变量表示将for循环随机抽取的区组组合以向量的形式存储,以便于后续输出结构化的随机分组表;
block_num变量用于记录区组的编号;
as.vector()函数是将block_set变量中存储的区组组合连接成与样本量相同的长度,以便与受试者编号对应;
data.frame()和write.csv()函数用于生成随机分组表并输出为csv格式,随机分组表中包含受试者编号、区组编号、分组。
实例二:为检验某一治疗肿瘤新化疗方案的临床疗效,选择常规化疗方案作为对照组,开展单中心双臂临床试验,需要受试者120名,每组各60名受试者,受试者入组期为3个月,化疗期为6个月,请按区组随机分组方法制定随机分组表。利用本文提供的自编函数block_random_m1()则可将参数size设定为120, 参数block_len设定为4,参数block_category设定为c(“T;T;C;C”, “T;C;T;C”, “T;C;C;T”, “C;C;T;T”, “C;T;C;T”, “C;T;T;C”),具体如下:
<-
res_block1 block_random_m1(size = 120,
block_len = 4,
block_category = c("T;T;C;C",
"T;C;T;C",
"T;C;C;T",
"C;C;T;T",
"C;T;C;T",
"C;T;T;C"))
函数输出的随机分组表见 表 tbl-random-2,其中共产生30个区组,每个区组内干预组和对照组各2名受试者。在实际的研究中区组长度的设定可适当变化,具体可见本文讨论部分。
kable(res_block1[c(1:12, 109:120), ])
ID | Block | Group | |
---|---|---|---|
1 | 1 | 1 | C |
2 | 2 | 1 | T |
3 | 3 | 1 | T |
4 | 4 | 1 | C |
5 | 5 | 2 | T |
6 | 6 | 2 | T |
7 | 7 | 2 | C |
8 | 8 | 2 | C |
9 | 9 | 3 | T |
10 | 10 | 3 | C |
11 | 11 | 3 | C |
12 | 12 | 3 | T |
109 | 109 | 28 | C |
110 | 110 | 28 | T |
111 | 111 | 28 | T |
112 | 112 | 28 | C |
113 | 113 | 29 | T |
114 | 114 | 29 | T |
115 | 115 | 29 | C |
116 | 116 | 29 | C |
117 | 117 | 30 | T |
118 | 118 | 30 | C |
119 | 119 | 30 | C |
120 | 120 | 30 | T |
6.3.2.2 方法二
第二种是先根据受试者数量生成随机数字,然后根据实验分组情况设定区组长度并依次对区组编号,再按区组编号分组对生成的随机数字编秩,最后根据各区组内随机数字秩的大小关系进行分组,其基本步骤为:1)按照实验设计的受试者样本量生成受试者编号;2)设定区组长度并根据受试者数量计算区组数并依次编号;3)利用随机函数生成与受试者编号数量相等的随机数字;4)按区组编号分组计算随机数字的秩;5)在各区组内,按照实验设计的组数及分组比例,对随机数字的秩序按照大小或者奇偶数进行分组;6)获得符合实验设计要求的随机分组表。本文同样提供了一段基于R语言的自编函数,能满足两组及多组且组间比例不同的区组随机分组需求,自编函数代码及具体含义如下:
library(dplyr, warn = FALSE)
<- function(size, block_len,
block_random_m2 grp = 2, T_2_C = "1:1"){
<- seq(1, size, 1)
id_num <- size / block_len
block
<- c()
block_num for (i in 1:block) {
<- append(block_num,
block_num rep(i, block_len))
}
<- runif(n = size,
random_seq min = 0, max = 1)
<- data.frame("ID" = id_num,
df "BlockNum" = block_num,
"RandomNum" = random_seq)
<- df %>%
df_rank group_by(BlockNum) %>%
mutate(Rank = rank(RandomNum))
<- as.numeric(substr(T_2_C, 1, 1))
ratio_T <- as.numeric(substr(T_2_C, 3, 3))
ratio_C
if (grp == 2){
$Group <-
df_rankifelse(df_rank$Rank <=
/(ratio_T + ratio_C),
block_len"T", "C")
else if (grp >2) {
} $Group <-
df_rankcut(df_rank$Rank,
breaks = grp,
labels = paste("Group", 1:grp))
}
# write.csv(df_rank,
# "block radomization table 2.csv",
# row.names = FALSE)
return(df_rank)
}
函数名为block_random_m2,函数共有4个参数,分别为size(样本量大小)、block_len(区组长度)、grp(分组数)和T_2_C(干预组与对照组比例),grp和T_2_C为默认参数,分别默认赋值为2和”1:1”;
seq()和runif()函数作用同上;
由于此方法需要进行分组编秩,考虑到R语言自带分组功能使用较为不便,本文在自编函数中引入了第三方程序包dplyr,可通过install.packages(“dplyr”)安装,其中group_by()、mutate()及%>%管道函数为dplyr包中引入的函数;
变量ratio_T和ratio_C,以及if条件语句中的功能与简单随机分组相同。
data.frame()和write.csv()函数用于生成随机分组表并输出为csv格式,随机分组表中包含受试者编号、区组编号、随机数字、分组。
同样以实例二为例,利用block_random_m2()函数可以获得同表2相似的结果,具体使用方式如下,输出结果见 表 tbl-random-3 。
<-
res_block2 block_random_m2(size = 120,
block_len = 4,
grp = 2,
T_2_C = "1:1")
kable(res_block2[c(1:12, 109:120), ])
ID | BlockNum | RandomNum | Rank | Group |
---|---|---|---|---|
1 | 1 | 0.2355342 | 3 | C |
2 | 1 | 0.1950892 | 1 | T |
3 | 1 | 0.7458419 | 4 | C |
4 | 1 | 0.2197261 | 2 | T |
5 | 2 | 0.0764369 | 1 | T |
6 | 2 | 0.6039926 | 4 | C |
7 | 2 | 0.4317742 | 2 | T |
8 | 2 | 0.4814425 | 3 | C |
9 | 3 | 0.7359029 | 2 | T |
10 | 3 | 0.8289157 | 3 | C |
11 | 3 | 0.9952405 | 4 | C |
12 | 3 | 0.5653650 | 1 | T |
109 | 28 | 0.7804675 | 4 | C |
110 | 28 | 0.6807008 | 3 | C |
111 | 28 | 0.3534151 | 2 | T |
112 | 28 | 0.3243371 | 1 | T |
113 | 29 | 0.3971316 | 2 | T |
114 | 29 | 0.5407747 | 4 | C |
115 | 29 | 0.4979685 | 3 | C |
116 | 29 | 0.1560650 | 1 | T |
117 | 30 | 0.6870941 | 1 | T |
118 | 30 | 0.7521457 | 2 | T |
119 | 30 | 0.8746505 | 4 | C |
120 | 30 | 0.7950460 | 3 | C |
6.3.3 分层随机分组(stratified randomization)
从操作层面而言,分层随机分组的本质是根据混杂因素对受试者进行分层后,再在各个层内采用简单随机或者区组随机的方法将受试者分配至不同的实验组。如在进行多中心肿瘤疾病的临床试验中,研究中心的情况和受试者肿瘤的分型分期会对治疗结局产生明显的影响,但是若采取简单随机或者区组随机分组的方法,很难确保不同研究中心入组的受试者特征相似,以及不同实验组间受试者肿瘤分型分期均衡,因此需要对此类明显可能造成实验结局偏倚的变量进行分层控制。假设存在3个研究中心,肿瘤共分为2型2期,那么就可以设置12个层,再在12个层内分别采用简单随机或者区组随机分组方法生成随机分组表,由于所采用方法同上,此处不再详细叙述。
6.4 讨论
本文详细梳理了当前阶段实验研究中常用的随机分组方法,并以此为基础,结合现阶段研究领域常用的开源统计软件R语言,详细介绍了实验研究中进行随机分组的具体过程,并通过自编函数的形式提出了集成化解决方案,且通过实例验证了自编函数的可靠性,可为研究者在开展实验研究的随机分组时提供借鉴和参考。有几点需要讨论说明,如下:
随机种子数的设定及管理。在RCT研究中,随机化通常和盲法联合使用,而随机分组的结果和规则是盲法的重要内容之一,因此除了需要保证不将随机分组结果透露至研究人员或者受试者外,还需要保证随机分组规则不可轻易被猜出。在利用计算机平台进行随机数生成时,合理设定随机种子数是保证生成的随机数可再现的重要前提。那么随机种子数的设定就十分关键,一般情况下随机种子数不可设置过于简单,如1或者2,亦或者年份。随机种子数应由项目组的统计人员保管,在实验揭盲或者结束之前不可透露给其他研究相关人员。
区组随机中区组长度的设定。区组长度由研究者确定,对于实验组为两组时通常取4和6(Kang, Ragan, and Park 2008)。考虑到区组长度较短且单一时,研究人员容易总结发现其中规律,可能会破坏盲法,因此在设定区组长度时可以采用复合方案,即在同一份随机分组表中同时设置4和6的区组长度,从而增加区组的复杂性,减少可预见性。另一方面,由于在具体的实验研究中,受试者脱落或者中途退出的情况并不少见,因此在生成随机分组表时,一般情况下推荐多设置1至2个区组,以备进行补充。
除本文中详细介绍的三种常用随机分组方法外,为了应对受试者入组的复杂情况,还有协变量自适应随机分组(covariate adaptive randomization)(Kang, Ragan, and Park 2008)和响应自适应随机分组( response-adaptive randomization)(刘晓燕 et al. 2008)方法,由于以上方法已有学者开发出了供R语言直接使用的软件包(carat),因此本文中未详细介绍,具体可见 <https://CRAN.R-project.org/package=carat>。