主页

原始数据 wget https://raw.githubusercontent.com/pzweuj/practice/master/R/Hist_PCA_Deno/Allen_BrainAtlas_12regions_Microarray.txt 柱状图 这里一共是有12个区域，画12个图并且放到一张图上。 brain12 <- read.table("Allen_BrainAtlas_12regions_Microarray.txt", header=TRUE) # 设置12种颜色 # 应该可以利用随机函数还有colors()函数来随机生成或者直接使用前12种颜色，更方便 col <- c("red", "blue", "or...

Plink是一款免费开源的用来分析全基因组的工具。 目前，已经更新了plink2。 PLINK的重点是对基因型/表型数据的分析，因此不支持此前的步骤（例如，研究设计和计划，从原始数据生成基因型或CNV调用）。 通过与gPLINK和Haploview集成，可以为后续的可视化，注释和结果存储提供支持。 我个人认为用新不用旧，因此这里会用plink2。 plink2可以在这里下载。windows的版本是编译好的exe文件。 操作简单。 具体的操作命令，可参阅General usage。 突出一个水字，因为时间不够用啦。

水一篇维持更新！ 最近几天发现我的WD服务变成了小红叉，而且点启动还会报错。在微软官网上看到的解决方案无法生效。 同时，我使用电脑更新的时候发现，报了0x800705b4错误。 赶紧去搜了一下，找到以下解决办法。 1 按Windows键+R，键入 services.msc 并回车，找到Windows Update服务，右击停止 2 在C:\Windows\SoftwareDistribution文件夹中，找到并删除Download和DataStore文件夹中的所有文件 3 然后重启Windows Update服务。 做完这些，我重启了一下电脑，然后更新不报错了，同时WD也自己启动了。 找一天有空再好好更新一篇。。。

背景 数据来源是精神病学基因组学联盟 可在这里申请下载数据。 参考文献PMID25056061 使用qqman包进行曼哈顿图等可视化分析。 安装R包 qqman install.packages("qqman") 载入与初始化 library("qqman") # 读入 SCZ <- read.table('rall.txt', header=TRUE, sep="\t") # 由于曼哈顿图需要的CHR变量是数字(int)，所以需要先把原始文件里的hg19chrc这一列改造一下 # 如把chr1改成1，把chrX改成23 # 去掉chr SCZ$hg19chrc <- gsub("chr", "", SCZ$hg19chrc) # 将X改成23 SCZ$hg1...

一般的，做完rna-seq的比对部分，就需要找出每个基因或者转录本的counts/RPKM/FPKM值。 之前介绍过featureCounts，统计counts非常快。 这次要说的是cufflinks，用来统计RPKM和FPKM值。 首先下载cufflinks wget http://cole-trapnell-lab.github.io/cufflinks/assets/downloads/cufflinks-2.2.1.Linux_x86_64.tar.gz tar -zxvf cufflinks-2.2.1.Linux_x86_64.tar.gz 然后把路径加入到环境中。 由于cufflinks需要的bam文件必须是排序过的，所以在采取hisat2进行比对的流程后，必须用s...

强推这篇文章，非常详尽，各种图例，眼花缭乱。 Analysis of single cell RNA-seq data 单细胞测序就是获取单个细胞遗传信息的测序技术，能够检出混杂样品测序所无法得到的异质性信息。 安装R包 SC3是一个非常好的用来分析单细胞测序的R包。 source("https://bioconductor.org/biocLite.R") biocLite("SC3") biocLite("scater") biocLite("SingleCellExperiment") 该R包有些参数是依赖perl的，所以还需要有perl的环境。 前数据处理 单细胞测序，如果原始数据只有一个fastq，要先把每个细胞的测序结果分开。 首先使用trim_galore来去接...

用clusterProfiler已经能出很好的图了。但是有时需要满足别的需求。 下载模拟的基因列表，也可以用自己的。 wget https://raw.githubusercontent.com/pzweuj/practice/master/R/GO_KEGG/gene.list 然后下面以GO CC注释为例 library(org.Hs.eg.db) library(clusterProfiler) library(ggplot2) # 读入文件 data <- read.table("gene.list",header=TRUE) data$GeneName <- as.character(data$GeneName) transID = bitr(data$G...

这个算是对上次ChIP-seq流程的一个补充。 在注释完peak之后，还可以进行找motif的操作。 蛋白质分子中的一些二级结构单元，往往有规则地聚集在一起形成全由α-螺旋、全由β-片层或α-螺旋与β-片层混合、均有的超二级结构基本形式，具体说，形成相对稳定的αα、βββ、βαβ、β2α和αTα等超二级结构又称模体(motif)或模序。几个motif可以组成功能单位(结构域)；一个或几个结构域构成活性中心(功能域)。 使用homer来找motif。 首先是安装 mkdir homer && cd homer wget http://homer.ucsd.edu/homer/configureHomer.pl perl configureHomer.pl -ins...