学术咖

Mindblown: a blog about philosophy.

  • 绘制盒形图

    输入表达矩阵数据文件 library(ggplot2) library(reshape2) exprSet <- read.table(file = “COS_deseq_counts_normalized.txt”, sep = “\t”, header = T) group_list <- factor(c(rep(“COS0”,3),rep(“COS1”,3),rep(“COS3”,3),rep(“COS6”,3),rep(“COS12”,3)), levels = c(“COS0”, “COS1″,”COS3″,”COS6″,”COS12″)) exprSet_L=melt(exprSet) colnames(exprSet_L)=c(‘Gene_ID’,’sample’,’value’) exprSet_L$group=rep(group_list,each=nrow(exprSet)) #以sample作为x轴，group作为填充颜色 ggplot(exprSet,aes(x=sample,y=value,fill=group)) + geom_boxplot() #以group作为x轴，sample作为填充颜色 #ggplot(exprSet,aes(x=group,y=value,fill=sample)) + geom_boxplot() # 修改sample名称 exprSet_L$sample=paste(exprSet_L$group,exprSet_L$sample,sep=”-“) ggplot(exprSet_L,aes(x=sample,y=value,fill=group)) + geom_boxplot() # 离群点的大小 cvbox=ggplot(exprSet_L,aes(x=sample,y=value,fill=group))+geom_boxplot(outlier.size=0) #outlier.size=-1 ，点将彻底消失 # 修改x轴字符方向 cvbox+theme(axis.text.x=element_text(angle=90,vjust=0.5,hjust=1)) # 修改图例位置和主题 cvbox+theme_bw()+theme(axis.text.x=element_text(angle=90,vjust=0.5,hjust=1),legend.position=”top”) # geom_jitter()：为图形设置抖动，防止相同点重叠 cvbox+theme_bw()+theme(axis.text.x=element_text(angle=90,vjust=0.5,hjust=1),legend.position=”top”) + geom_jitter()

  • 叔本华：作为意志和表象的世界（一）

    第一篇  世界作为表象初论 　服从充分根据律的表象 经验和科学的客体 跳出童年时代吧，朋友，觉醒呵！ ——J·J·卢梭 “世界是我的表象”：这是一个真理，是对于任何一个生活着和认识着的生物都有效的真理；不过只有人能够将它纳入反省的，抽象的意识罢了。并且，要是人真的这样做了，那么，在他那儿就出现了哲学的思考。于是，他就会清楚而确切地明白，他不认识什么太阳，什么地球，而永远只是眼睛，是眼睛看见太阳；永远只是手，是手感触着地球；就会明白围绕着他的这世界只是作为表象而存在着的，也就是说这世界的存在完全只是就它对一个其他事物的，一个进行“表象者”的关系来说的。这个进行“表象者”就是人自己。如果有一真理可以先验地说将出来，那就是这一真理了，因为这真理就是一切可能的、可想得到的经验所同具的那一形式的陈述。它比一切，比时间、空间、因果性等更为普遍，因为所有这些都要以这一真理为前提。我们既已把这些形式都认作根据律的一些特殊构成形态，如果其中每一形式只是对一特殊类型的表象有效，那么，与此相反，客体和主体的分立则是所有那些类型的共同形式。客体主体分立是这样一个形式：任何一个表象，不论是哪一种，抽象的或直观的，纯粹的或经验的，都只有在这一共同形式下，根本才有可能，才可想象。因此，再没有一个比这更确切，更不依赖其他真理，更不需要一个证明的真理了；即是说：对于“认识”而存在着的一切，也就是全世界，都只是同主体相关联着的客体，直观者的直观；一句话，都只是表象。当然，这里所说的对于现在，也对于任何过去，任何将来，对于最远的和近的都有效；因为这里所说的对于时间和空间本身就有效；而又只有在时间、空间中，所有这些[过去、现在、未来、远和近]才能区别出来。一切一切，凡已属于和能属于这世界的一切，都无可避免地带有以主体为条件[的性质]，并且也仅仅只是为主体而存在。世界即是表象。 这个真理决不新颖。它已包含在笛卡儿所从出发的怀疑论观点中。不过贝克莱是断然把它说出来的第一人；尽管他那哲学的其余部分站不住脚，在这一点上，他却为哲学作出了不朽的贡献。康德首先一个缺点就是对这一命题的忽略，这在本书附录中将有详尽的交代。与此相反，吠檀多哲学被认为是毗耶舍的作品，这里所谈的基本原理在那里就已作为根本命题出现了，因此印度智者们很早就认识这一真理了。威廉·琼斯在他最近《论亚洲哲学》（《亚洲研究》，第四卷第164页）一文中为此作了证，他说：“吠檀多学派的基本教义不在于否认物质的存在，不在否认它的坚实性、不可入性、广延的形状（否认这些，将意味着疯狂），而是在于纠正世俗对于物质的观念，在于主张物质没有独立于心的知觉以外的本质，主张存在和可知觉性是可以互相换用的术语。”这些话已充分地表出了经验的实在性和先验的观念性两者的共存。 在这第一篇里，我们只从上述的这一方面，即仅仅是作为表象的一面来考察这世界。至于这一考察，虽无损于其为真理，究竟是片面的，从而也是由于某种任意的抽象作用引出来的，它宣告了每一个人内心的矛盾，他带着这一矛盾去假定这世界只是他的表象，另一方面他又再也不能摆脱这一假定。不过这一考察的片面性就会从下一篇得到补充，由另一真理得到补充。这一真理，可不如我们这里所从出发的那一个，是那么直接明确的，而是只有通过更深入的探讨，更艰难的抽象和“别异综同”的功夫才能达到的。它必然是很严肃的，对于每一个人纵不是可怕的，也必然是要加以郑重考虑的。这另一真理就是每人，他自己也能说并且必须说的：“世界是我的意志。” 在作这个补充之前，也就是在这第一篇里，我们必须坚定不移地考察世界的这一面，即我们所从出发的一面，“可知性”的一面：因此，也必须毫无抵触心情地将当前现成的客体，甚至自己的身体（我们就要进一步谈到这点）都仅仅作为表象看，并且也仅仅称之为表象。我们希望往后每一个人都会确切明白我们在这样做的时候，只仅仅是撇开了意志；而意志就是单独构成世界另外那一面的东西；因为这世界的一面自始至终是表象，正如另一面自始至终是意志。至于说有一种实在，并不是这两者中的任何一个方面，而是一个自在的客体（康德的“自在之物”可惜也不知不觉的蜕化为这样的客体），那是梦呓中的怪物；而承认这种怪物就会是哲学里引人误入迷途的鬼火。 转载：叔本华哲学智慧

  • GSEA命令行

    前边已经说到MSigDB数据库的entrez ID和symble ID都是人源的，其他物种要想做GSEA的话，就必须要ID转换，做成属于自己研究物种的数据库，这一步是关键。想要根据简单的Nr注释或其他注释来转换是行不通的，比如我研究的物种为硬骨鱼类，Nr注释结果分值最高的前20个一般也都是鱼类，就算注释到了人源的基因，也可能存在版本不对的问题。因此我们要知道MSigDB数据库基因ID的版本，然后下载对应的全基因组蛋白序列，然后比对卡阈值，做成自己的数据库。 GSEA的运行可以界面操作和命令行操作，界面操作教程很多就不多做赘述，它的不足之处除了操作繁琐之外，里面几个较大的库需要的内存较高，像我的4G渣渣电脑，好几个都运行失败，不得不转到服务器上运行；而且界面也不能批量运行。 java -cp gsea-3.0.jar -Xmx10000m xtools.gsea.Gsea -res ma/normalized_counts_all.gct -cls ma/samples.cls#COS1_versus_COS0 -gmx ma/ma.h.all.v6.1.symbols.gmt -collapse false -mode Max_probe -norm meandiv -nperm 1000 -permute gene_set -rnd_type no_balance -scoring_scheme weighted -rpt_label COS1vsCOS0_h -metric Ratio_of_Classes -sort real -order descending -create_gcts false -create_svgs false -include_only_symbols true -make_sets true -median false -num 100 -plot_top_x 20 -rnd_seed timestamp -save_rnd_lists false…

  • GSEA分析流程

    什么是GSEA？ 基因集富集分析 (Gene Set Enrichment Analysis, GSEA)用来评估一个预先定义的基因集的基因在与表型相关度排序的基因列表中的分布趋势，从而判断其对表型的贡献。基因集合富集分析检测基因集合而不是单个基因的表达变化，因此可以包含这些细微的表达变化，预期得到更为理想的结果。 GSEA原理 给定一个排序的基因表L和一个预先定义的基因集S (比如编码某个代谢通路的产物的基因, 基因组上物理位置相近的基因，或同一GO注释下的基因)，GSEA的目的是判断S里面的成员s在L里面是随机分布还是主要聚集在L的顶部或底部。这些基因排序的依据是其在不同表型状态下的表达差异，若研究的基因集S的成员显著聚集在L的顶部或底部，则说明此基因集成员对表型的差异有贡献，也是我们关注的基因集。 GSEA分析中几个关键的概念 计算富集得分 (ES, enrichment score)，ES反应基因集成员s在排序列表L的两端富集的程度。计算方式是，从基因集L的第一个基因开始，计算一个累计统计值。当遇到一个落在s里面的基因，则增加统计值。遇到一个不在s里面的基因，则降低统计值。每一步统计值增加或减少的幅度与基因的表达变化程度（更严格的是与基因和表型的关联度）是相关的。富集得分ES最后定义为最大的峰值。正值ES表示基因集在列表的顶部富集，负值ES表示基因集在列表的底部富集。 评估富集得分(ES)的显著性，通过基于表型而不改变基因之间关系的排列检验 (permutation test)计算观察到的富集得分(ES)出现的可能性。若样品量少，也可基于基因集做排列检验 (permutation test)，计算p-value。 多重假设检验矫正，首先对每个基因子集s计算得到的ES根据基因集的大小进行标准化得到Normalized Enrichment Score (NES)。随后针对NES计算假阳性率。（计算NES也有另外一种方法，是计算出的ES除以排列检验得到的所有ES的平均值） Leading-edge subset，对富集得分贡献最大的基因成员。 MSigDB数据库 MSigDB（Molecular Signatures Database）数据库中定义了已知的基因集合：http://software.broadinstitute.org/gsea/msigdb 包括H和C1-C7八个系列（Collection），每个系列内容为：H: hallmark gene sets（效应）特征基因集合，共50组（比如细胞凋亡特征基因集）；C1: positional gene sets 位置基因集合，根据染色体位置，共326个；C2: curated gene sets：（专家）共识基因集合，基于通路、文献等（这部分包括我们熟悉的KEGG信号通路）；C3: motif gene sets：模式基因集合，主要包括microRNA和转录因子靶基因两部分；C4: computational gene sets：计算基因集合，通过挖掘癌症相关芯片数据定义的基因集合；C5: GO gene sets：Gene Ontology 基因本体论，包括BP（生物学过程biological process，细胞原件cellular component和分子功能molecular…

  • DESeq分析流程

    source(“https://bioconductor.org/biocLite.R”) biocLite(“DESeq”) library(“DESeq”) database <- read.table(file = “macrophage_genes.count_table.matrix”, sep = “\t”, header = T, row.names = 1) countData <- database[,1:3] condition <- factor(c(“A”,”B”,”C”))#无生物学重复 #type <- factor(c(rep(“A”,3), rep(“B”,3))) 有生物学重复 database <- round(as.matrix(countData))#取整数型 dds <- newCountDataSet(database,condition) dds <- estimateSizeFactors(dds) dds <- estimateDispersions(dds, method=”blind”, sharingMode=”fit-only” )#无生物学重复 dds <- estimateDispersions(dds) # 有生物学重复 resAvsB <- nbinomTest(dds,”A”,”B”) resAvsC <- nbinomTest(dds,”A”,”C”) table(resAvsB$pval…

  • Python replace()方法

    描述 Python replace() 方法把字符串中的 old（旧字符串） 替换成 new(新字符串)，如果指定第三个参数max，则替换不超过 max 次。 语法 replace()方法语法： str.replace(old, new[, max]) 参数 old — 将被替换的子字符串。 new — 新字符串，用于替换old子字符串。 max — 可选字符串, 替换不超过 max 次 返回值 返回字符串中的 old（旧字符串） 替换成 new(新字符串)后生成的新字符串，如果指定第三个参数max，则替换不超过 max 次。 实例 以下实例展示了replace()函数的使用方法： #!/usr/bin/python str = “this is string example….wow!!! this is really string”; print str.replace(“is”, “was”); print str.replace(“is”, “was”, 3); 以上实例输出结果如下：…

  • ggplot2：气泡图及散点图小结

    气泡图也属于散点图的一种，在散点图的基础上改变点的形状，大小和颜色 1.如何改变点的形状 用自带的mtcars演示 p = ggplot(mtcars,aes(wt,mpg)) p + geom_point(shape=17,size=4) p + geom_point(aes(shape = factor(cyl)),size=4) 2.改变点的大小：气泡图 加载包、数据的读入及数据格式 library(ggplot2) pathway = read.table(“R0-vs-R3.path.richFactor.head20.tsv”,header=T,sep=”\t”) 3.初始化数据 pp = ggplot(pathway,aes(richFactor,Pathway)) pp + geom_point() 气泡图（三维数据）,size 可以是一个值，也可以是数据中的一列 pp + geom_point(aes(size=R0vsR3)) 添加颜色变化——四维数据 scale_colour_gradient：自定义一个连续型的配色，常用于热图； pbubble = pp + geom_point(aes(size=R0vsR3,color=-1*log10(Qvalue))) pbubble pbubble + scale_colour_gradient(low=”green”,high=”red”) expression函数改变样式，[]是用来添加下标，^是用来添加上标 pr = pbubble + scale_colour_gradient(low=”green”,high=”red”) + labs(color=expression(-log[10](Qvalue)),size=”Gene number”,x=”Rich factor”,y=”Pathway name”,title=”Top20 of…

  • ggplot2：火山图

    火山图是散点图的一种，也是描点，用 geom_point()绘图 1.加载包，读入数据及数据格式 library(ggplot2) data = read.table(“R0-vsR3.isoforms.filter.tsv”,header=T,row.names=1) 2.绘图 r03 = ggplot(data,aes(log2FC,-1*log10(FDR))) r03 + geom_point() 3.改变点的颜色 r03 + geom_point(color=”red”) r03 + geom_point(aes(color=”red”)) r03 + geom_point(aes(color=significant)) # 按照“ significant”这一列定义点的颜色； 4.设置坐标轴范围和标题 函数xlim()，ylim()，规定X、 Y轴范围；labs(title=“..”,x=“..”,y=“..”)，确定标题和坐标轴标签；expression函数改变样式，[]是用来添加下标，^是用来添加上标 r03xy = r03 + geom_point(aes(color=significant)) + xlim(-4,4) + ylim(0,30) r03xy + labs(title=”Volcano plot”,x=”log2(FC)”) r03xy + labs(title=”Volcano plot”,x=expression(log[2](FC)),y=expression(-log[10](FDR))) 5.自定义颜色 scale_color_manual(): 自定义颜色配色； r03xyp = r03xy + labs(title=”Volcano…

  • ggplot2：散点图

    1.加载ggplot2包，读取数据及数据格式 library(ggplot2) fpkm = read.table(“all.fpkm”,header=T,row.names=1) head(fpkm,10) 2.绘图 mp = ggplot(fpkm,aes(C1_FPKM,C2_FPKM)) mp + geom_point() 3.取log10后重新作图 mp = ggplot(fpkm,aes(log10(C1_FPKM),log10(C2_FPKM))) mp + geom_point() 4.添加直线 slope设置斜率，intercept设置截距，color设置线条颜色，size设置线条粗细 mp + geom_abline(slope=1,intercept=0,color=”red”,size=2) + geom_point() ggplot2的核心思想就是图层，可以像PS一样一层层的添加图层，以加号连接，更详细参数可参考ggplot2的帮助文档http://docs.ggplot2.org/current/ 视频教程：http://www.omicshare.com/class/home/index/singlev?id=34

  • WGCNA分析流程

    WGCNA是一个R包，对一个完全不会R的人来说，确实费了一番功夫，不过也将我对R的学习提前提上日程。 分析步骤： 1.WGCNA安装 source(“http://bioconductor.org/biocLite.R”) biocLite(c(“AnnotationDbi”, “impute”, “GO.db”, “preprocessCore”)) install.packages(“WGCNA”) 2.输入数据准备 原始数据https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE48213，共56个样本；对数据的合并可参考http://www.biotrainee.com:8080/thread-603-1-1.html，输入文件格式如下：   setwd(‘WGCNA/’) fpkm=read.table(‘test.txt’,sep = ‘\t’,stringsAsFactors = F) head(fpkm) dim(fpkm) ##共56个细胞系，36953个基因的表达矩阵 names(fpkm)##细胞系名称 3.数据预处理 library(reshape2) library(‘WGCNA’) RNAseq_voom <- fpkm WGCNA_matrix = t(RNAseq_voom[order(apply(RNAseq_voom,1,mad), decreasing = T)[1:3000],]) ##top 3000 mad genes，t参数将矩阵变为符合WGCNA要求的形式：行名为gene，列名为样品 datExpr <- WGCNA_matrix save(datExpr, file = “AS-green-FPKM-01-dataInput.RData”) 4.确定最佳beta值   powers = c(c(1:10), seq(from = 12, to=20, by=2))…

Got any book recommendations?

Get In Touch