par(mar = c(4, 4, 0.5, 0.5), mfrow = c(1, 2))
set.seed(123)
plot(pval <- t(replicate(1000, {
x1 = rnorm(100, mean = 0, sd = runif(1, 0.5, 1))
x2 = rnorm(100, mean = 1, sd = runif(1, 2, 5))
c(t.test(x1, x2, var.equal = TRUE)$p.value, t.test(x1, x2,
var.equal = FALSE)$p.value)
})), xlab = "P-value: equal variance", ylab = "P-value: unequal variance",
pch = 20, asp = 1)
abline(0, 1)
plot(pval[, 1], pval[, 2] - pval[, 1], xlab = "P-value: equal variance",
ylab = "Diff of p-values (unequal var - equal var)", pch = 20)
过程是：从两个正态总体中生成样本，第一个总体均值为0，标准差随机取自U(0.5, 1)，第二个总体均值为1，标准差取自U(2, 5)，显然两个总体标准差不相等，那么在t检验时设定和不设定方差相等的选项对结果有多大影响？把两种情况的P值都画出来：左图是原始P值，可见基本在对角线上，说明大致相等，若眼神儿不好，可看右图，即P值的差异，可见方差不等时P值偏大（原因很简单，因为Welch校正的自由度小于等于不校正的自由度，样本量相等的时候统计量的分母即标准误一样，因此统计量完全一样，自由度越小，P值越大嘛），但大多少呢？其实也没大多少。
方差齐与不齐时t检验的结果对照
近来有些奇怪，有几位R core们居然给本小子写邮件，让本小子着实感到吃惊。比如，首先是Brian Ripley，这是R core中的core，前面提到过他在R源代码中的突出贡献，这位真人不露面、网上找不着照片的大佬，前段时间给我发封邮件说，你小子的animation包的启动消息不正规啊，因为我用suppressPackageStartupMessages()无法屏蔽启动消息；我一看，R里面居然还存在这么长名字的函数，顺便学习了message()函数，从此不再用老土的cat()函数了，后来考虑了一下，干脆把启动消息去掉了，library(animation)不会再有任何提示消息。
然后是我发现Duncan Temple Lang这位不靠谱的大叔做着一些我很喜欢的不靠谱的事情，于是乎对Omegahat心向往之，一来二去聊了聊，将来有机会一定要会一会他。
Sweave对注释的处理是要么完全去除，然后R代码会被整理整齐，要么完全保留，但R代码也保留原样，而我一直希望既能保留注释又能整理代码，这才诞生了animation包中tidy.source()函数（在小邱聪明的技巧下），前段时间想想给Friedrich Leisch，也就是Sweave的作者，发封邮件说了这个事情，打探一下是否能多设置一些Sweave选项，比如把parse()和deparse()函数以选项的形式抽象出来，这样就可以实现既整理代码又保留注释的功能了，不过大叔貌似很忙，回了一封邮件就再也没有音信了，后来由于Michael Friendly对Sweave的一些功能请求在R-help上发了邮件，我们一干人等通过Duncan Murdoch间接了解到Friedrich的确很忙，不过好消息是圣诞节过后Sweave可能会有更新，届时用户可以自行设置图形设备，不必局限在PDF和EPS。但整理代码的事情仍然遥遥无期……唉，还得用硬性Hack的方法。
Martin Maechler前面提过，看到我们开R会，说要向The R Journal交报告啊，回头再跟他谈谈明年R会议的事情。
最意想不到的是，Duncan Murdoch刚才居然给我发个邮件问问题，额滴神啊，这位大叔可是Rtools的管理者、若干个包的作者（rgl等）啊。不过大叔问的是Flash的问题，还好我知道那么一点点，算是能解决。趁此机会，干脆回问两个C语言问题，子曾经曰过：问一个够本，问两个赚一个。
这极北苦寒之农村今年比较反常，据说往年都是感恩节一定会下雪，而今年就没下。刚纳闷儿咋12月还不下雪的时候，便下开了。不过雪也不大，地上只是铺了薄薄一层，比今年北京那场雪差远了。昨夜回家路上，想起小学时咿咿呀呀背的：
日暮苍山远，天寒白屋贫。
柴门闻犬吠，风雪夜归人。
第二届中国R语言会议的北京会场过几个小时就要开幕了。这次会议比上次的准备更加匆忙，大约也就只有一个多月时间准备，但大家都很卖力，在此先感谢一下各位组织者：邱怡轩、张翔、焦静、陈堰平、范建、蒋安华以及关菁菁；说起小邱同学，如我上次所说，我真是有点怕给这位拼命三郎安排任务，从别人口中了解到他为这次会议每天马不停蹄焦头烂额四处奔波，我心中甚为感叹；张翔呢，我没想到他会担起这次会议组织者的角色，上海会场在他的带领下也办得有声有色（看看会议通知页面的宣传海报多么亮丽），和焦静两人拉赞助、发传单、安排吃住，作为已经工作的人，对一门自由软件如此费心，甚为难得；焦静呢，现在不在统计专业（生态），却帮忙做着一门统计软件的推广，跑校区、定会场、找领导，忙得不亦乐乎；陈堰平作为R的老用户挑起大梁，相信经过上次植物所培训一战，对这次会议的组织应该更有把握；fan版主也是位拼命三郎，COS论坛招生就业版自他上任之后所有帖子和资料被整理得井井有条，使你不得不敬佩，这年头能如此发狠的人不多见啊；关菁菁同学嘛，说实话刚开始小邱介绍的时候我想不起来她是哪一级了（为啥我总觉得她是研究僧呢），上一届R会议她参加了，而且中午没去吃饭，留在会场帮我们看东西，这次又主动提出愿意帮忙组织，我自然是很高兴。对于参会者诸如魏太云以及刘思喆和李舰二位大师兄的献计献策一并致谢。
这次会议有不少去年的熟面孔（如丁鹏、左辰、王化儒、奚潭等），新参加的人里面有我认识和不认识的，报名演讲的名单也给了我很多惊喜，比如钟其顶，算是一个老朋友了，三年前在我一次R报告的时候就认识了，后来我们一直用R做一些食品行业的应用，效果还是很不错的，尤其是今年初几位师弟在我的牵线搭桥下过去实习之后，挖掘了更多R的应用价值；再比如那场“地质环境调查监测研究中的R应用”（作者来自中国地质环境监测院，政府机构下属事业单位），其图形着实让我惊叹了一番，没想到R在这样的单位已经被人研究到了这种程度（R和Google Earth都用上了），太出乎我的意料了，看来我的统计图形书可以放到更开阔的边界上把各种稀奇古怪的应用都介绍一下；再比如陈丽云，这位以技术派面目出现的lady，要来讲讲计量，想当年，被本小子一句玩笑代码惹得好奇心起，装了R，然后被打击了一番；再比如孙晓燕，最后关头杀了过来，不知道是不是被李晓煦老师给“忽悠”的；还有中科院的WebR，相信也是很价值的应用；上海会场请到了汤银才老师，这位也是中国R语言的一位元老人物，想想四五年前网上一搜就是汤老师的那个PPT。
一个月前我往R-help发了个会议通知，前几天又补充了一下会议内容，R core之一Martin Maechler注意到我们的会议，给我发了封邮件说希望这次会议能写一篇报告发给The R Journal，正中下怀，本小子就是这么打算的。我想明年第三届R会议就放在暑假开好了，到时候请一些R core成员过来讲讲课什么的，应该也是很有可能的。
前两天给吴老写了封邮件说起这事，吴老曰：
我是一个行将退出战场的老兵，我想说的是：
祝贺第二届中国R语言大会胜利召开！
开放、绿色、功能强大、具有源源不断巨大资源的R不仅有必要而且一定能够在中国推广和发展。
吴老是第一位把R引进人大统计学院的老师，此后他的弟子们也纷纷用R，这才有了我接触R的机会。
昨天Ripley教授向R提交了第50000次修改，Romain Francios对SVN的日志数据做了一些简单分析，我个人一直关心Ripley是不是整天不用睡觉（你看这老爷子一天到晚都在邮件列表中出没），这次正好验证一下他是不是24小时工作，数据和R代码参见50000 Revisions Committed to R。
Ripley从1999年加入R核心团队，从上图可以看出，他显然是不需要睡觉的——每个小时都可能有commit。时间分布呈双峰：早上7点到10点、下午3点到6点。看看其他人的工作时间，很容易发现Martin倾向于早起干活，而Peter倾向于每天晚上12点之后干活。
Romain的博客中有SVN的日志数据可以下载，感兴趣的同志们可以继续分析R core的工作行为。
上周上课，我们老爷子又说SAS is extremely powerful，SAS很靠谱，就差明说SAS没有Bug了。我正在整理课程笔记，SAS的事情，我改天要去找老爷子好好谈一谈。R从不说自己能担保什么，大家拼命找bug，拼命改进，这是开源软件的共同特征——没有人付钱，但就是有一群疯子半夜3点还在写代码。商业软件一向说自己能保证什么，可是一个bug二十年都没人会去修正（例），还好我们看不见源代码，要是能看见，后果不堪设想。
前些日子有位童鞋在R-help邮件列表里问如何生成那种单词大小与其出现频率成比例的图，这玩意儿也就是通常所说的标签云（Tag Cloud）。我琢磨了一下WordPress的插件wp-cumulus，发现其原理很简单，不过就是将标签信息以XML形式通过JavaScript传递给一个Flash文件，所以也很容易用R去实现这个传递过程，即：将文本、超级链接以及频数写成XML，然后嵌入到HTML文件中。整个过程参见Creating Tag Cloud Using R and Flash / JavaScript (SWFObject)这篇日志，函数源代码和示例数据都可以从那里下载。
效果是这样的：
（通过RSS阅读的童鞋们请打开原文链接在浏览器中观看，否则啥都看不到）
x区间上取一系列点，然后计算y = f(x)，把这些点描在笛卡尔坐标系中，然后用铅笔手工连起来。这个办法很笨很原始，但很有效，奇怪的是，到计算机如此发达的今天，仍然有人不断问如何用软件画曲线，仿佛中学时学的东西都还给数学老师了。例如，若干天前，有人问我如何用R画这个曲线：
我只好提醒他回忆中学的描点法。R里面的gamma()和curve()函数是现成的，把数学公式写成相应的代码就可以了：
curve((sqrt(2) * gamma(x/2) * 2)/(sqrt(x - 1) * gamma(x - 1)), 2, 50) # 当然，如果你非要按照x = seq(2, 50, ...)然后计算f(x)然后plot(x, f(x), type = "l") # 那也可以，只不过curve()就是干这事的，不必麻烦分两步走了
描点法画函数曲线
如同某些网文说人一辈子的守则在幼儿园都已经教完了，就看你长大还记不记得。所以每当面对一个貌似复杂的问题时，首先要想想幼儿园有没有学过。
“武功再高，也怕菜刀”绝对是真理，哈哈。
最近到处看到关于第一期The R Journal上线的消息，随便举几个例子，比如魏太云那嘎达、刘思喆那嘎达、陈钢那嘎达、Zhanwu Dai的邮件、Paulo那嘎达（blogspot又一次被和谐了）、David Smith那嘎达等等。
这次R Journal的面世，在R界是一件重大事情，R News这个刊名确实有点小家山寨，而其内容已经堪比正规统计刊物，所以更名The R Journal是理所当然的。之所以大家全都提到了这件事情我还要写这篇博文是因为：【插播广告：有哪位同学/老师愿意承担整理论文集的任务，请速速与我联系。这个任务包括：催促作者们完善论文（比如我个人的论文至今还没写完）、联系出版社、召集评委、排版等。谢谢！】
- 我们第一届中国R语言会议的会议纪要发表在第一期R Journal上了，这是很有历史意义的。尽管当初主编大人说i suspect this will be a nice article，后来我也没收到其它邮件，于是应验了no news is good news，发表了。这里要向各位客官交代的是，我们的论文集在我的电脑里扔了半年也没整理，出版社也没再联系，抱歉抱歉。
- 我在这一期The R Journal上看到一个人：Xuefei Mi。遥想当年我在多特蒙德参加useR! 2008的时候，某一天下午散会之后在多特蒙德大学外面等火车，看见一位童鞋，看长相我琢磨着可能是日本人吧，他也瞅瞅我，最终开口了，原来都是中国人……在火车上简单聊了聊，得知他和Torsten Hothorn有某种关系（导师？忘了），如今在R Journal上又看见这位兄台的名字了。这世界是不是很小呢？
- 前两条都还忍住了，当遇到第三件事情的时候，终于决定操起键盘写一下了。以前我曾经提到过Felix Andrews，也就是
playwith和latticist包的作者，与Deepayan Sarkar合作了latticeExtra包。Felix今天给我写了个邮件，算是回复了我174天前给他的GTalk留言（邀请他来参加R会议），说在R Journal上看见会议纪要了，当时很忙所以没去参加，一堆英文下来，最后附了一句中文：秀月也好了。她每天全天说话，说得很流利。她喜欢我们家里的猫。他家小女儿还是那么可爱，哈哈，刚到Felix的网站找了找，可惜没找着照片（只有她在她妈妈肚子里时的照片）；Felix的中文也很流利了。话说Felix回到ANU建了一个堪培拉的R Group，我看他们每个月都组织活动，其实我们也可以组队，每次说这个事情都是没有带头人……眼看我也要撤了，希望北京能赶快诞生一个定期聚会的R小组啊。
总结：R在中国的繁荣昌盛不会太远了，但同志们需要大大努力。
此时此刻，俺哆嗦着告诉各位客官，俺真的中奖了。
话说两个月前，我介绍了一下John Chambers软件奖，当时并没有指望获奖，只是打探打探情况，搞明白之后等我明年好好琢磨一个靠谱的软件之后再正式参加。结果，不知是美国人民都忙着对付金融危机了，还是一共就我一个人参加，或是怎么的，今日打开邮箱，发现我被通知得奖（中奖？）了。我愣了半天，觉得这是在做梦。
好吧，两个月前说了要请客，那么那篇日志底下的观众朋友们请赶紧拿出小本本记下，“谢某人欠俺一顿大餐”。
做动画都能获奖……敢问谱在何方啊……管它呢，反正是美国人的钱，俺遵旨，8月去华盛顿领钱，顺便领奖。
话说前天在海关开会的时候没事干，便偷偷在下面写去年R会议的报告，三下五除二写完昨天发给The R Journal了，今天主编大人说，i suspect this will be a nice article，看来悬了，估计是发不出来了。R core咋就不能照顾照顾广大中国R人民呢……
最近评论