R软件 蒙特卡罗模拟
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
R使用指南
打开R
下图是R软件的主窗口,R软件的界面与Windows的其他编程软件类似,由一些菜单和快捷按钮组成。快捷按钮下面的窗口便是命令输入窗口,它也是部分运算结果的输出窗口,有些运算结果则会在新建的窗口中输出。
当一个R 程序需要你输入命令时,它会显示命令提示符。默认的提示符是>。技术上来说,R 是一种语法非常简单的表达式语言(expression language)。它大小写敏感,因此A 和a 是不同的符号且指向不同的变量。可以在R 环境下使用的命名字符集依赖于R 所运行的系统和国家(就是系统的locale 设置)。通常,数字,字母,. 和都是允许的(在一些国家还包括重音字母)。不过,一个命名必须以. 或者字母开头,并且以. 开头时第二个字符不允许是数字。基本命令要么是表达式(expressions)要么就是赋值(assignments)。如果一条命令是表达式,那么它将会被解析(evaluate),并将结果显示在屏幕上,同时清空该命令所占内存。赋值同样会解析表达式并且把值传给变量但结果不会自动显示在屏幕上。命令可以被(;)隔开,或者另起一行。基本命令可以通过大括弧(f和g) 放在一起构成一个复合表达式(compound expression)。注释几乎可以放在任何地方7。一行中,从井号(#)开始到句子收尾之间的语句就是注释。如果一条命令在一行结束的时候在语法上还不完整,R 会给出一个不同的提示符,默认是+。该提示符会出现在第二行和随后的行中,它持续等待输入直到一条命令在语法上是完整的。该提示符可以被用户修改。在后面的文档中,我们常常省略延续提示符(continuation prompt),以简单的缩进表示这种延续。
R的帮助
首先先来看看如何使用帮助文件
这里有两个方式:
在R中输入
help.start()
通过启动HTML形式的在线帮助(使用你的计算机里面可用的浏览器)。你可以用鼠标点
击上面的链接。
如图:
英文足够好的话,这就是R的最佳教材(可惜我不行!!!)。
当你看到某个不知道语句可以通过这方式来找到答案:
例如plot()语句,你想知道他是干什么的那么就在R中输入“?plot()”
看到了吗,很详细的解释。
看懂了plot()是做什么的了吗?,没看懂看下面的例子。
下面的这部分会话让你在操作中对R 环境的一些特性有个简单的了解。你对系统的许多特性开始时可能有点不熟悉和困惑,但这些迷惑会很快消失的。
登录,启动你的桌面系统。
R 程序开始,并且有一段引导语。
(在R 里面,左边的提示符将不会被显示防止混淆。)
help.start()
启动HTML 形式的在线帮助(使用你的计算机里面可用的浏览器)。你可以用鼠标点
击上面的链接。
最小化帮助窗口,进入下一部分。
x <- rnorm(50)
y <- rnorm(x)
这里x<-rnorm(50)产生50 个标准正态数据,y<-rnorm(x)是产生和x为数一样多的一组标准正态数据。
plot(x, y)
画二维散点图。一个图形窗口会自动出现。
ls()
查看当前工作空间里面的R 对象。
rm(x, y)
去掉不再需要的对象。(清空)。
x <- 1:20
等价于x = (1; 2; : : : ; 20)。
w <- 1 + sqrt(x)/2
标准差的`权重'向量。
dummy <- data.frame(x=x, y= x + rnorm(x)*w)
dummy创建一个由x 和y构成的双列数据框,查看它们。fm <- lm(y ~ x, data=dummy)
summary(fm)
拟合y 对x 的简单线性回归,查看分析结果。
fm1 <- lm(y ~ x, data=dummy, weight=1/w^2)
summary(fm1)
现在我们已经知道标准差,做一个加权回归。
attach(dummy)
让数据框中的列项可以像一般的变量那样使用。
lrf <- lowess(x, y)
做一个非参局部回归。
lines(x, lrf$y)
增加局部回归曲线。
abline(0, 1, lty=3)
真正的回归曲线:(截距0,斜率1)。
编写自己的函数
R 语言允许用户创建自己的函数(function)对象。R 有一些内部函数可以用在其他的表达式中。通过这个过程,R 在程序的功能性,便利性和优美性上得到了扩展。学写这些有用的函数是一个人轻松地创造性地使用R 的最主要的方式。
需要强调的是,大多数函数都作为R 系统的一部分而提供,如mean(),
var(),postscript()等等。这些函数都是用R 写的,因此在本质上和用户写的没有差别。一个函数是通过下面的语句形式定义的,
> name <- function(arg 1 , arg 2 , ...) expression
其中expression是一个R 表达式(常常是一个成组表达式),它利用参数argi计算最终的结果。该表达式的值就是函数的返回值。
可以在任何地方以name(expr1 , expr2 , ...) 的形式调用函数。
计算机模拟方法Monte Carlo方法
下面主要介绍最基本的计算机模拟方法Monte Carlo方法,此方法的基本思想是将各种随机时间的概率特征(概率分布、数学期望)与随机事件的模拟联系起来,用试验的方法确定事件的相应概率与数学期望.因而, Monte Carlo方法的突出特点是概率模型的解是由试验得到的,而不是计算出来的。此外,模拟任何一个实际过程, Monte Carlo方法都需要用到大量的随机数,计算量很大,人工计算是不可能的,只能在计算机上实现.
在概率论中,著名的Buffon掷针问题就是用统计试验的方法求圆周率π的典型代表.现在用Monte Carlo方法模拟一下Buffon投针试验,对于该问题大家都已经了解了其解法,下面就不再赘述,我们只从随机模拟的角度度该问题进行求解。
针与平行线相交的充分必要条件是
sin
2
l
xθ
≤
。Buffon的投针试验在计算机上实
现,需要一下两个步骤:
产生随机数.首先产生n个相互独立的随机变量θ,x的抽样序列
,,1,2,, i i
x i n θ=,
其中
(0,),(0,)
2
i i
a
U x U
θπ
.
模拟检验:检验不等式
sin
2
i i
l
xθ
≤
是否成立.
若上式成立,表示第i次试验成功(即针与平行线相交).设n次试验中有k次成
功,则π的估值为
2
ˆl n
ak
π
*
=
,其中a l>,均为预先给定。
将上述步骤编写成R模拟程序。模拟代码如下:
buffon<-function(n, l=0.8, a=1) #定义一个新的函数buffon,包括三个自变量(n ,l, a)
{