“统计之都”站的“网站导读”栏目的设立是为了以简短的形式向大家介绍一些有意思而且有水平的统计学文章，不求理论之复杂，但求统计学之生活化，让大家看到一些统计学的“另类”面目。若这个目的达不到，那么我希望大家读完这些导读文章之后能说一句“哇，原来统计不是会计啊/不是做报表的啊/不是数学啊”也足够了。另外，现在网上很多文章都是抄来抄去（更恶劣的是不加出处的抄袭），我们觉得这种做法极其无聊，是对原始作者的极大不尊重，也容易造成以讹传讹误导不明真相的围观群众，本站这个栏目的建立，也是基于这一点考虑之上提供一种“引用他人文章”的示例，很傻很天真地希望互联网的抄袭现象能够有所收敛。言归正传： 一、用数字作弊要小心：关于内贾德的选票 日前伊朗的大选可谓轰轰隆隆，颇引人注目，内贾德在胜出之后却引来一片质疑，今日又有传闻说穆萨维才是真正的胜出者。总之疑云重重，那么让政治家玩政治家的游戏吧，我们从另一个视角来关注一下这次选举。密歇根大学的教授Walter R. Mebane, Jr.这几天一直在分析选票数据，今天的文章参见：Note on the presidential election in Iran, June 2009（注意论文和数据以及R代码都在更新中，如果不能访问，请到他的主页上找） 文章主要基于两点理论去检验选票数据： 本福特定律（Benford’s Law）：生活中的数据里，1~9这9个数字在首位的出现并非均匀分布，例如1出现在一个数字的首位的概率约为1/3，而不是想象中的1/9，越往后的数字出现概率越低。不过，作者得到的选票数据是按地区汇总的，而对于汇总数据，我们往往难以发现它作弊的嫌疑，因为汇总数据倾向于符合本福特定律，从这一点上，作者没有找到足够的证据证明选票数据作弊； 检验离群点：手头的数据只有选票数，没有其它变量，这种情况下统计建模的局限性很大，似乎各种工具都施展不开，不过作者还是在“艰苦”的条件下建立了过度散布的二项回归模型（overdispersed binomial regression），因变量是二分变量（内贾德 vs 穆萨维），自变量是选票数，看看哪些地区的选票数有离群点出现，熟悉伊朗政局情况的人从这些离群点可能发现违背常理的现象（选票过高？过低？） 看来统计学家总是有办法拷问数据，不知咱国内是否有这样的统计学家呢？ 想起本科上课的时候我们有一位老师提到陈毅元帅的一句话：“莫伸手，伸手必被捉。”在数据上，还是莫造假吧，只要有人较真，一定有办法看出来的。 二、飞机被流星撞到的概率多大：法航失事之后的计算 前几天，本站作者、COS论坛元老刘思喆用R计算了一下飞机被流星击中而失事的概率，发现一架飞机在11小时飞行过程中被流星击中的概率是双色球中一等奖概率的1/100，看样子可以舒口气，不过坏消息是，20年中被击中的概率就陡然上升到5%了。看来还是双腿走路可靠…… 然而我想知道的是： 过去20年有飞机被流星击中过么？ Poisson分布合适么？或者，流星砸向地球的分布是什么？有没有天文学相关的实证呢？ 三、让R帮你查看是否中了500万大奖：“懒惰”的彩民 这一篇依旧是把统计和R融入生活的刘思喆：500万？去买双色球！（被评论为标题党）彩民们每天眼巴巴对着彩票网以及自己的彩票的若干位数字看，实在是很辛苦，因此思喆老大体谅大家，挑灯夜战写了一段R代码，按照该文的描述，彩民们每天只要打开R，就知道自己有没有中奖。我提议，要是哪位读者因此中奖了一定要给思喆老大提成以及给统计之都捐赠！ 我后来的想法是，干脆写个批处理文件（R CMD BATCH），开机自动运行好了，连R都不用打开。懒惰是创新之源，这话没错。 R还有很多“奇怪”的应用，在作者的主页贝吉塔行星中能找到更多。 四、其它文章 画函数曲线本来是一件简单的事情，为什么我们总是忘记初中的“描点法”呢：画曲线的通用办法：描点法画图（谢益辉） 标签云将词语文本的大小与其某种属性（例如重要性、出现频率）关联起来，因此标签云图可以直观展示一些词语的属性。例如文本越大表示出现频率越高，那么一眼看去，最大的词就是最热频词了，这里介绍了一种用R生成标签云的方法：Creating Tag Cloud Using R and Flash / JavaScript (SWFObject)（谢益辉） 最后强烈建议对统计软件感兴趣（尤其是对R）的同志们订阅Journal of Statistical Software的RSS，你会经常发现一些稀奇古怪的论文：http://www.jstatsoft.org/rss（啥？你到现在还不知道RSS是神马东西？快来人把这个火星观众撵出去） 若您平时看到觉得写得好的文章，请勿忘推荐给我们，联系邮箱：contact@cos.name。谢谢！
阅读全文 »

在研究两个变量的关系时，一般会先看看它们的散点图，在图中两变量的关系还是比较直观的，大致可以判断是否线性相关及相关性大小如何，是否是非线性相关等。而到底什么是相关呢？相关其实就是知道一件事对了解另一件事的帮助的大小。实际中，如果对某一事物不太了解，但是对与其有一定联系的另一事物有所了解，如果这种联系很强，那我们对于那件不了解的事物就有了更多的信息，或者说对这个不了解的事物有了更大的自信去预测。其实这也是研究中的一种常用的方法。
阅读全文 »

新华网刚发布了一个关于学生冬季长跑的调查结果（于2009年4月27日13:52访问），一共调查了100人，结果中却出现了92.79%这样的比例数字，有常识的读者都知道，世上不存在0.79个人，因此这里面必然有某个地方是错的（姑且不妄言造假）。这则消息让我马上想起《统计陷阱》这本书，我们生活中有多少陷阱呢？ 从消息编辑人员的角度来说，他们可能觉得保留两位小数显得“精确”，而这种“精确精神”从数学的角度来说显得既可爱又可笑，如果小数位能表示精确，那何不保留100位小数呢？ 从统计人的角度来看，对这种调查报告中的比例数字应该有足够的警觉。很多调查报告并不会告诉我们究竟样本量多大（在这一点上新华网的调查网还比较诚实），这种情况下，我们应该弄清究竟调查了多少对象，当样本量很小的时候我们会怀疑这个调查的代表性。当我们看到比例66%的时候也许能想起来这是2/3（猜测样本量是3的倍数），但对29.1667%这个比例我们未必能很快反应出来分子和分母是多少，若报告公布方没有说明样本量，我们只能自己猜测；对于667这样的数字，我们很容易猜测这是6循环的四舍五入。最终大概思路就是拿比例去挨个乘以一系列整数，看看哪个结果接近整数，从而“还原”原来的分式n/N。以下是简单的R代码： > digit = ((1:100) * 0.29166666)%%1 # 整除1之后的“余数” > plot(digit, ylim = c(0, 1)) > idx = which((1 - digit) < 1e-05 | (digit - 0) < 1e-05) # 与0或1很靠近时 > points(idx, digit, pch = 20) > abline(v = idx, lty = 2) > axis(3, idx) > idx * 0.29166666 7 14 21 28 我们很容易发现分母（样本量）是24的倍数，因为上图中24的倍数乘以29.167%得到的结果很靠近整数；而具备某种特征的样本数量为7的倍数。根据具体的调查背景，我们可以自己猜测报告方究竟调查了多少人：24人？太少；960人？为什么不是1000人？…… 上面只是统计侦查的小游戏而已，当我们具备更多统计知识储备之后，便可以去考虑一些具体的统计模型输出是否存在造假嫌疑。我想，P值在0.05之下且很靠近0.05的时候，或P值一律接近于0的时候，我们不妨以小人之心揣测这个模型也许有问题。当因子分析中，50个变量能根据载荷被准确划分到作者预先设定的5个因子中时（5列因子载荷一律都是只在某个因子上取值极大），这个分析也许存在嫌疑。当然，所有的“小人之心”的前提假设都是：理想情况在现实中是不容易出现的（这是赤裸裸的假设检验逻辑）。
阅读全文 »

大家好，为了促进大家对统计之都的了解，并锻炼各位会员的统计应用能力，即日起我们推出“COS竞赛”系列活动。第一期活动的主要任务是分析统计之都英文网站（http://cos.name/en/）的会员数据，从中找出识别正规会员和机器人（垃圾、广告、自动注册）会员的规律。 数据背景 原始数据来自phpBB论坛的phpbb_users数据库，其中包含用户id、用户名、是否激活、Email、发帖数等字段，其中我们要研究的因变量是“是否激活”（user_active），它取值0和1，分别代表该用户是否被管理员激活，1表示是正规会员，0表示该会员是机器人或不知如何激活的人类。 我们对原始数据首先做了隐私处理，删掉了一些敏感字段，如用户名等，然后对Email和MSN帐号做了如下处理： x = read.csv("cos_en_users.csv", stringsAsFactors = FALSE, encoding = "UTF-8") x$user_email_nchar = nchar(x$user_email) x$user_email = sub("^.*@", "", x$user_email) x$user_msnm_nchar = nchar(x$user_msnm, allowNA = TRUE) x$user_msnm = sub("^.*@", "", x$user_msnm) write.csv(x, "cos_en_spam.csv", row.names = FALSE) 即：去掉了帐号中@字符前面的字符串，只剩下域名字符串，然后在原数据中添加了两个变量*_nchar表示Email和MSN帐号的字符数。 数据下载 整理之后的数据下载：统计之都英文网站会员数据； 所有变量名的解释参见：phpbb_users的数据表结构说明；注意其中所有时间数据都是Unix时间戳格式，起点为1970-01-01 00:00:00，可以用R函数as.POSIXlt(, origin = "1970-01-01 00:00:00")转化为日期。 分析目的 找出有效的规则区分会员类型：是否正规会员。本数据是一个高度不平衡的数据，取值为1的会员非常少： > table(x$user_active) 0 1 25911 92 这给判别分析带来了不小的难度，如：即使将所有正规会员判别为机器人，正确率也在99.6462%，因此仅仅看模型的正确率可能不是合适的评判标准。我们最终考察的指标包括：预测正确率、模型简洁性、程序效率、分析结果新颖性。 注意事项 使用数据时请注意数据有效性，可以从原始信息中生成新的变量用来做预测： user_session_page：该变量的取值意义不太确定，建议不使用； user_lastvisit：变量取值不稳定，建议不使用，user_session_time可能更适合于作为用户最后的访问时间； 由于非激活用户不能发帖，因此他们的user_posts必然都是0；而我们的分析目标是在用户注册之后马上就能从注册信息获知是否机器人注册，所以建议建模时不要使用这个变量； user_style取值只能为1，是个常数，因此不必使用该变量； 建议着重分析用户注册信息中的签名档、邮箱、网站链接，可以从这些文本数据生成新的变量，如邮箱域名是否以”.ru”（俄罗斯）结尾，签名档是否含有”free”等具有垃圾特征的词汇，等等； 参加方式 本次竞赛以邮箱投稿的方式接收作品，请将您的作品发送至contact cos.name（替换为@），作品应满足以下要求： 注明使用软件工具的详细信息（版本、操作系统等）； 评委可在其它地方相同条件下重复您的完整分析过程，包括：数据转换、建模、预测、图表的输出；为了保证预测结果的可重复性，请适当增大预测中估计正确率的交叉验证次数（或其它验证方法的次数）； 写明最终的分析结论； 文档格式不限，但有如下优先考虑顺序：Sweave > LaTeX/LyX...
阅读全文 »

用GeoBUGS做简单的空间数据分析 第一节 实例介绍基本的空间模型 GeoBUGS是WinBUGS的一个模块，专门用来分析空间数据（spatial data)。由于和WinBUGS的基本模型结合得比较好，所以被广泛地使用。目前的GeoBUGS除了自身的地图格式外，还支持Splus, ArcInfo 以及 EpiMap的地图格式。当然了，在使用的时候需要做适当的转化才行。 下面是一个简单的例子，大家也可以在GeoBUGS的Manual中找到它。模型假设为条件自回归模型 Conditional Autoregressive（CAR）。数据为苏格兰唇癌疾病数据，反映的是苏格兰56个郡的唇癌发病率。这个数据比较经典，Clayton and Kaldor (1987) 和 Breslow and Clayton (1993)都曾在他们的论著中分析过该数据。
阅读全文 »

第一节 WinBUGS数据分析案例 在这一节中，我将拿一个经典的研究数据，利用WinBUGS给出简单的分析。首先介绍一下这个数据：Seeds seed O. aegyptiaco 75 seed O. aegyptiaco 73 Bean Cucumber Bean Cucumber r n r/n r n r/n r n r/n r n r/n 10 39 0.26 5 6 0.83 8 16 0.5 3 12 0.25 23 62 0.37 53 74 0.72 10 30 0.33 22 41 0.54 23 81 0.28 55 72 0.76 8 28 0.29 15...
阅读全文 »