R语言-随机数和抽样
前言
在Matlab、R或者S-PLUS等软件中做随机数模拟时,经过会遇到set.seed()这个函数。随机数的产生需要有一个随机的种子,因为用计算机产生的随机数是通过递推的方法得来的,必须有一个初始值。用同一台电脑,且在初始值和递推方法相同的情况下,可以产生相同的随机序列。
用计算机产生的是“伪随机数”。用投色子计数的方法产生真正的随机数 , 但电脑若也这样做 , 将会占用大量内存;用噪声发生器或放射性物质也可产生真正的随机数,但不可重复。而用数学方法产生最适合计算机,这就是周期有限,易重复的 “伪随机数”。
set.seed函数括号内的参数是一个数,可能会经常看到该数取不同的值。这个函数的目的主要是为了产生相同的随机数,这样一来,每次运行的结果就是一样的,这使得结果可重复被验证。参数的选择很随意,如2、23、124等等。
为了验证上面的结论,我分别取 set.seed(2)、 set.seed(2)、 set.seed(3),生成一组随机变量,对生成的随机变量作图,由此验证了前两组随机变量在种子数相同的情况下,得到的值是一致的(黑、红),而第三组情况(蓝)则与前两者不同。
目录
1. 随机数生成
2. 抽样模拟
3. 设置随机数种子
1. 随机数生成
作为一种统计分析语言,R是一个生成各种统计分布功能随机数的综合性图书馆。在这篇文章中,我想专注于这个简单的问题:我如何生成一个随机数?
答案取决于你想要什么样的随机数生成?让我们通过例子说明。
a) 生成 5.0 和 7.5 之间的随机数
如果你想生成一个十进制数规定的最低和最高之间的任何值(包括分数值)同样是可能的,使用runif()功能。 这个函数生成均匀分布的值。 这里是如何生成一个 5.0 和 7.5 之间的随机数的方法:
> x1 <- runif(1, 5.0, 7.5)
> x1
[1] 5.573882
当然,当你运行这个,你会得到一个不同的数字,但它一定会在 5.0 和 7.5 之间。 你不会得到准确值 5.0 或 7.5 。
如果你想生成多个随机的值,不要使用一个循环。您可以生成多个值一次通过指定您要作为第一个参数runif()值的数目。这里是如何产生10个 5.0 和 7.5 之间的值:
> x2 <- runif(10, 5.0, 7.5)
> x2
[1] 6.871893 6.518833 5.372662 6.808643 5.490282 7.381824
[7] 7.476985 5.569176 6.182591 5.285149
b) 生成 1 到 10 之间的随机整数
这看起来像最后一个相同的运动,但现在我们只希望完整的数字,而不是分数值。 为此,我们使用的示例函数:
> x3 <- sample(1:10, 1)
> x3
[1] 3
c) 各种分布的随机数生存函数
rnorm(n, mean=0, sd=1) # 正态分布
rexp(n, rate=1) # 指数
rgamma(n, shape, rate=1, scale=1/rate) # r 分布
rpois(n, lambda) # 泊松
rt(n, df, ncp) # t 分布
rf(n, df1, df2, ncp) # f 分布
rchisq(n, df, ncp=0) # 卡方分布
rbinom(n, size, prob) # 二项分布
rweibull(n, shape, scale=1) # weibull 分布
rbata(n, shape1, shape2) # bata 分布
2. 抽样模拟
a) 在 1 到 10 之间随机抽取 1 个数
> x1 <- sample(1:10, 1)
> x1
[1] 3
第一个参数是一个有效的数字向量生成(这里的数字1到10),第二个参数表示应返回一个数字。
b) 在 “一组数” 之间随机抽取 多 个数
如果我们要生成多个随机数,我们必须增加一个额外的参数,表示允许重复:
# 有放回抽取
> x2 <- sample(1:10, 5, replace=T)
> x2
[1] 1 5 9 9 8
如果我们要生成多个随机数,我们必须增加一个额外的参数,表示不允许重复:
# 无放回抽取
> x3 <- sample(1:40, 6, replace=F)
> x3
[1] 3 40 20 19 28 11
c) 对一组数据进行乱序排序(随机抽取,顺序随机)
你可以使用同样的想法产生的任何载体的随机子集,甚至不包含数字。 例如,选择10个不同的美国各州随机:
> sample(state.name, 10)
[1] "Delaware" "Vermont" "Rhode Island" "Tennessee"
[5] "Arizona" "Mississippi" "Virginia" "Alaska"
[9] "Georgia" "Louisiana"
> sample(state.name, 52)
Error in sample.int(length(x), size, replace, prob) :
cannot take a sample larger than the population when 'replace = FALSE'
> length(state.name)
[1] 50
> sample(state.name, 50)
[1] "Oregon" "Georgia" "Maine"
[4] "Idaho" "Alaska" "Tennessee"
[7] "Indiana" "Wyoming" "Montana"
[10] "Utah" "Florida" "North Carolina"
[13] "Nevada" "Virginia" "Pennsylvania"
[16] "North Dakota" "Wisconsin" "Alabama"
[19] "Mississippi" "New Hampshire" "Delaware"
[22] "Arizona" "Massachusetts" "Vermont"
[25] "Maryland" "Missouri" "Michigan"
[28] "Connecticut" "Colorado" "New York"
[31] "Oklahoma" "California" "Washington"
[34] "South Carolina" "West Virginia" "Hawaii"
[37] "Illinois" "Arkansas" "Kentucky"
[40] "Iowa" "Kansas" "Rhode Island"
[43] "New Mexico" "South Dakota" "New Jersey"
[46] "Nebraska" "Ohio" "Texas"
[49] "Louisiana" "Minnesota"
> sample(state.name, 50)
[1] "Arkansas" "California" "Texas"
[4] "Tennessee" "Montana" "Massachusetts"
[7] "North Dakota" "Oregon" "Delaware"
[10] "Hawaii" "South Dakota" "Connecticut"
[13] "South Carolina" "Kansas" "Washington"
[16] "West Virginia" "Georgia" "Maine"
[19] "Wyoming" "Illinois" "Nebraska"
[22] "Idaho" "Maryland" "Ohio"
[25] "Indiana" "Louisiana" "Mississippi"
[28] "Iowa" "New York" "Minnesota"
[31] "Rhode Island" "New Mexico" "New Jersey"
[34] "Pennsylvania" "Kentucky" "Utah"
[37] "Nevada" "Arizona" "Missouri"
[40] "North Carolina" "New Hampshire" "Wisconsin"
[43] "Alaska" "Vermont" "Alabama"
[46] "Florida" "Oklahoma" "Michigan"
[49] "Colorado" "Virginia"
>
d) 随机从矩阵(数据框)中选取一部分对象
-----------------------------------------------------
(col.name=colnames(mtcars))
(row.name=rownames(mtcars))
#列名向量不返回抽样
(sam.col.name=sample(col.name,10,replace=FALSE))
#行名向量不返回抽样
(sam.row.name=sample(row.name,10,replace=FALSE))
B=mtcars[sam.row.name,sam.col.name]
B
-----------------------------------------------------
> (col.name=colnames(mtcars))
[1] "mpg" "cyl" "disp" "hp" "drat" "wt" "qsec" "vs" "am"
[10] "gear" "carb"
> (row.name=rownames(mtcars))
[1] "Mazda RX4" "Mazda RX4 Wag"
[3] "Datsun 710" "Hornet 4 Drive"
[5] "Hornet Sportabout" "Valiant"
[7] "Duster 360" "Merc 240D"
[9] "Merc 230" "Merc 280"
[11] "Merc 280C" "Merc 450SE"
[13] "Merc 450SL" "Merc 450SLC"
[15] "Cadillac Fleetwood" "Lincoln Continental"
[17] "Chrysler Imperial" "Fiat 128"
[19] "Honda Civic" "Toyota Corolla"
[21] "Toyota Corona" "Dodge Challenger"
[23] "AMC Javelin" "Camaro Z28"
[25] "Pontiac Firebird" "Fiat X1-9"
[27] "Porsche 914-2" "Lotus Europa"
[29] "Ford Pantera L" "Ferrari Dino"
[31] "Maserati Bora" "Volvo 142E"
> #列名向量不返回抽样
> (sam.col.name=sample(col.name,10,replace=FALSE))
[1] "mpg" "qsec" "wt" "disp" "carb" "cyl" "vs" "am" "hp"
[10] "gear"
> #行名向量不返回抽样
> (sam.row.name=sample(row.name,10,replace=FALSE))
[1] "Merc 230" "Pontiac Firebird" "Maserati Bora"
[4] "Hornet Sportabout" "Merc 450SE" "Merc 280"
[7] "Duster 360" "Merc 280C" "Dodge Challenger"
[10] "Camaro Z28"
> B=mtcars[sam.row.name,sam.col.name]
> B
mpg qsec wt disp carb cyl vs am hp gear
Merc 230 22.8 22.90 3.150 140.8 2 4 1 0 95 4
Pontiac Firebird 19.2 17.05 3.845 400.0 2 8 0 0 175 3
Maserati Bora 15.0 14.60 3.570 301.0 8 8 0 1 335 5
Hornet Sportabout 18.7 17.02 3.440 360.0 2 8 0 0 175 3
Merc 450SE 16.4 17.40 4.070 275.8 3 8 0 0 180 3
Merc 280 19.2 18.30 3.440 167.6 4 6 1 0 123 4
Duster 360 14.3 15.84 3.570 360.0 4 8 0 0 245 3
Merc 280C 17.8 18.90 3.440 167.6 4 6 1 0 123 4
Dodge Challenger 15.5 16.87 3.520 318.0 2 8 0 0 150 3
Camaro Z28 13.3 15.41 3.840 350.0 4 8 0 0 245 3
>
-----------------------------------------------------
e) 有放回的随机抽样示例
在e里随机有放回地取N数,形成一个向量
e <- c(1:6)
N <- 100
s <- sample(x=e, size=N, replace=TRUE)
table(s)
> e <- c(1:6)
> N <- 100
> s <- sample(x=e, size=N, replace=TRUE)
> table(s)
s
1 2 3 4 5 6
16 19 14 21 19 11
f) 用sample函数实现随机排序(抽样)
sample函数第一个参数x表示被抽样的向量,第二个参数size表示抽取的样本个数,第三个参数replace表示是不是有放回的抽样。
-----------------------------------------------------
Vec <- 1:10
# 无放回抽样
sample(x = Vec, size = 5)
sample(x = Vec, size = 5)
# 只要第二个参数size的数值和第一个参数x的长度相同就可以实现随机排序。
sample(x = Vec, size = 10)
sample(x = Vec, size = 10)
# 有放回抽样
sample(x = Vec, size = 10, replace = TRUE)
sample(x = Vec, size = 10, replace = TRUE)
-----------------------------------------------------
> Vec <- 1:10
> # 无放回抽样
> sample(x = Vec, size = 5)
[1] 9 3 1 5 4
> sample(x = Vec, size = 5)
[1] 3 1 4 9 10
> # 只要第二个参数size的数值和第一个参数x的长度相同就可以实现随机排序。
> sample(x = Vec, size = 10)
[1] 1 6 2 10 7 8 3 9 5 4
> sample(x = Vec, size = 10)
[1] 6 9 10 7 1 5 8 3 4 2
> # 有放回抽样
> sample(x = Vec, size = 10, replace = TRUE)
[1] 3 6 2 2 8 8 3 6 10 8
> sample(x = Vec, size = 10, replace = TRUE)
[1] 7 5 10 5 7 3 1 6 4 7
-----------------------------------------------------
g) rep()和replicate()批量取随机数
问题:假设我想从符合正态分布的数据集中随机抽取2个数据,然后排序,而这样的数据我需要20对,你会怎么做?
# 方法1:rep()函数
rep(sort(sample(rnorm(n=100,mean = 0,sd = 1),2)),20)
> rep(sort(sample(rnorm(n=100,mean = 0,sd = 1),2)),20)
[1] -0.7433281 0.1923704 -0.7433281 0.1923704
[5] -0.7433281 0.1923704 -0.7433281 0.1923704
[9] -0.7433281 0.1923704 -0.7433281 0.1923704
[13] -0.7433281 0.1923704 -0.7433281 0.1923704
[17] -0.7433281 0.1923704 -0.7433281 0.1923704
[21] -0.7433281 0.1923704 -0.7433281 0.1923704
[25] -0.7433281 0.1923704 -0.7433281 0.1923704
[29] -0.7433281 0.1923704 -0.7433281 0.1923704
[33] -0.7433281 0.1923704 -0.7433281 0.1923704
[37] -0.7433281 0.1923704 -0.7433281 0.1923704
# 方法2:replicate()函数
replicate(n=20,expr=sort(sample(rnorm(n=100,mean = 0,sd = 1),2)))
> replicate(n=20,expr=sort(sample(rnorm(n=100,mean = 0,sd = 1),2)))
[,1] [,2] [,3] [,4]
[1,] 0.1214492 0.1682693 -1.4858967 -0.2435413
[2,] 0.5112295 0.6430365 0.1062271 0.1608918
[,5] [,6] [,7] [,8]
[1,] -0.7942942 -0.9092859 -0.38338123 -1.3953220
[2,] 0.6556288 1.0844620 -0.06132436 -0.3722372
[,9] [,10] [,11] [,12]
[1,] -0.6628561 -0.9390822 -1.772896 -1.886010
[2,] 0.6303991 1.2778667 -1.588206 -1.099241
[,13] [,14] [,15] [,16]
[1,] -1.1779682 -0.6320272 -0.7623491 -0.03874195
[2,] -0.5253614 0.1441117 0.4199327 0.70510390
[,17] [,18] [,19] [,20]
[1,] 0.1193869 -1.774997 -1.4304266 -0.7024574
[2,] 0.3463737 -1.338911 0.1309915 0.3383929
>
小结:rep()返回的是向量,replicate()返回的是矩阵。
下面列出两个函数的用法:
rep():
rep(x, ...)
rep.int(x, times) # 每个元素重复次数
rep_len(x, length.out) # 生成向量长度
# 随机数生成器
replicate(),replicate(n, expr, simplify = "array")
3. 设置随机数种子
> set.seed(1234)
> x <- rnorm(10)
> x
[1] -1.2070657 0.2774292 1.0844412 -2.3456977 0.4291247
[6] 0.5060559 -0.5747400 -0.5466319 -0.5644520 -0.8900378
>
> # 设置了随机数种子后,每次产生的分布数都是相同的
> set.seed(1234)
> x <- rnorm(10)
> x
[1] -1.2070657 0.2774292 1.0844412 -2.3456977 0.4291247
[6] 0.5060559 -0.5747400 -0.5466319 -0.5644520 -0.8900378
>
> # 移除了随机数种子后,产生的随机分布改变了
> rm(.Random.seed)
>
> x <- rnorm(10)
> x
[1] 0.6977298 0.6580453 0.5422976 0.1453601 -0.4882267
[6] 0.6037503 -0.8052218 0.3463612 -2.3368599 0.6643264
>
©哈尔滨商业大学 银河统计工作室
银河统计工作室成员由在校统计、计算机部分师生和企业数据数据分析师组成,维护和开发银河统计网和银河统计博客(技术文档)。专注于数据挖掘技术研究和运用,探索统计学、应用数学和IT技术有机结合,尝试大数据条件下新型统计学教学模式。