本文心得自：The Split-Apply-Combine Strategy for Data Analysis, Hadley Wickham, Journal of Statistical Software, April 2011, V.40.
引子：
我们常常会遇到这样的问题，数据量很大，并不需要依顺序来依次处理。合理分块处理，并最终整合起来是一个不错的选择。这也就是所谓的Split-Apply-Combine Strategy策略。这在速度上会有比做一个loop有优势，因为它可以并行处理数据。
什么时候我们需要使用到化整为零的策略呢？有以下三种情况：
- 数据需要分组处理
- 数据需要按照每行或者每列来处理
- 数据需要分级处理，和分组很类似，但是分级时需要考虑分级之间的关系。
化整为零策略有点类似于由Google推广的map-reduce策略。当然map-reduce策略的基础是网格，而这里的Split-Apply-Combine的基础完全可以是单机，甚至不支持并行处理的单机都可以。
然而，化整为零并不是一个很直观的编程过程。最直观的过程是使用Loop循环。这里使用一个例子来讲解一下如何实现化整为零策略。在plyr包中有数据ozone,它是一个三维矩阵(24X24X72)，其中最后一维72是指的6年12个月每个月的结果。也就是ozone是一个包括了连续72个月24Ｘ24的三维矩阵数据。三维分别是lat,long,time。我们需要由对时间robust linear model之后的残基residuals。
> library(plyr) # need for dataset ozone > library(MASS) # need for function rlm > month < - ordered(rep(1:12, length=72)) #set time sequence > #try one set > one < - ozone[1,1,] > model < - rlm(one ~ month - 1); model Call: rlm(formula = one ~ month - 1) Converged in 9 iterations Coefficients: month1 month2 month3 month4 month5 month6 month7 month8 month9 month10 month11 month12 264.3964 259.2036 255.0000 252.0052 258.5089 265.3387 274.0000 276.6724 277.0000 285.0000 283.6036 273.1964 Degrees of freedom: 72 total; 60 residual Scale estimate: 4.45 > deseas < - resid(model)
UTF8_E[......]
Tags: bioconductor, 教程, 编程
引子
闲来无事，翻看CELL杂志，发现很多基因组的图都使用circos来作图，于是就去circos.ca网上看了一眼，发现果然是个基因组研究绘图强大工具，应该为生物信息员掌握。于是花了点时间，来试验它的每一个设置。上网搜索发现，它的中文资料少得可怜，于是将心得体会做一总结，形成这一系列教程。希望能对急于掌握circos又不擅长阅读英语的人有所帮助。－－糗世界之糗糗
下载与安装
下载地址：http://circos.ca/software/download/circos/
circos是基于perl的脚本程序。它的安装难度在于安装好perl以及它所需要的模块。对于windows用户，可以试着安装Strawberry Perl或者ActiveState Perl。这两者都是不错的选择。对于Unix/linux/MacOS用户，很可能你已经安装了perl，否则，你可以到http://www.perl.org/get.html去下载安装。
我们需要测试一下perl的环境， UNIX/Linux/MacOS用户
> which perl /usr/bin/perl > perl -v This is perl, v5.10.0 built for ...
Windows用户
> perl -v This is perl, v5.10.0 built for ...
接着，我们将下载下来的circos程序解压，假设它的目录是circos-x.xx
> cd circos-x.xx > bin/circos -man
你可能得到一个帮助页面，那么你安装circos已经成功。也可能得到是出错信息。[......]
首先声明，这里并不是我用perl来写个算法。如果需要算法相关的东西，得去看文献。
这里，其实只是借对格兰氏阳性细菌蛋白质定位的预测为例，来说明如何利用他人提供的在线查询功能来实现批处理。
假设现在有一个完整的基因组需要您去预测它当中的每一个开放阅读框翻译出来的蛋白质可能的细胞内定位，怎么办呢？手工一个一个提交到网站上去？一共会有四五千个蛋白，等你提交完，你的手和大脑都会不工作了吧？要不自己下载个软件来本地预测吧？我试过去安装那些要求的软件环境，也许是我的系统过新吧，一个c语言库的版本，一个g77让我就头大得不知道该怎么继续下去。于是我还是下定决心，用perl的lwp来伪装成浏览器提交申请，自动批处理吧。也许一觉醒来，全部都做完了。
我们要用到的网站是http://www.psort.org/psortb/。据网站上宣传，Based on a study last performed in 2010, PSORTb v3.0.2 is the most precise bacterial localization prediction tool available. 第二个原因就是can currently submit one or more Gram-positive or Gram-negative bacterial sequences or archaeal sequences in FASTA format。这对于研究格兰氏阳性菌的我来说的确很不错的网站。
最基本的，用lwp来虚拟提交一份表单上去：
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 | #!/usr/bin/perl use LWP::UserAgent; $ua = LWP::UserAgent->new; $ua->agent("Mozilla 8.0 beta"); use HTTP::Request::Common qw(POST); my $req = (POST 'http://www.psort.org/psortb/results.pl', ["seq" => ">$genes{name}\n$genes{translation}", "organism" => "bacteria", "gram"=>"positive", "advancedgram"=>"none", "format"=>"long", "sendresults"=>"display"); $request = $ua->request($req); print $request->as_string; |
结果得到的，怎么都是500 Internal Server Error。于是在网上狂google，也没有找到直接的答案。半夜两点，突然想起，为什么不自己想法来解决问题。于是开始用tcpdump抓包，对比从firefox发送出去的包和perl发送出去的包，具体看看有什么不同。[......]
近期评论