基因数据处理77之从vcf文件中提取某条染色体的数据

1.代码:

/**
  * @author xubo
  */
package org.gcdss.cli.vcf

import org.apache.spark.{SparkConf, SparkContext}

/**
  * Created by xubo on 2016/5/23.
  */
object extractGRCH38chr20vcf {
  def main(args: Array[String]) {
        val conf = new SparkConf().setAppName(this.getClass().getSimpleName().filter(!_.equals('$'))).setMaster("local[4]")
//    val conf = new SparkConf().setAppName(this.getClass().getSimpleName().filter(!_.equals('$')))
    val sc = new SparkContext(conf)

    println("start:")
    val startTime = System.currentTimeMillis()
           val rdd = sc.textFile(args(0)).filter { each =>
      val line = each.split("\t")
      if (line(0) == args(1).toString) {
        true
      } else {
        false
      }
    }

    println(rdd.count())

    val loadTime = System.currentTimeMillis()
    println("load time:" + (loadTime - startTime) + " ms")

    //    rdd.repartition(1).saveAsTextFile(args(0) + "chr" + args(1) + ".vcf")

    val saveTime = System.currentTimeMillis()
    println("save time:" + (saveTime - loadTime) + " ms")
    println("run time:" + (saveTime - startTime) + " ms")
    sc.stop
  }
}

2.脚本:

hadoop@Master:~/xubo/project/vcf/extractGRCH38chr20vcf$ cat run.sh 
    #!/usr/bin/env bash  
    spark-submit   \
--class  org.gcdss.cli.vcf.extractGRCH38chr20vcf \
--master spark://219.219.220.149:7077 \
--conf spark.serializer=org.apache.spark.serializer.KryoSerializer \
--conf spark.kryo.registrator=org.bdgenomics.adam.serialization.ADAMKryoRegistrator \
--jars /home/hadoop/cloud/adam/lib/adam-apis_2.10-0.18.3-SNAPSHOT.jar,/home/hadoop/cloud/adam/lib/adam-cli_2.10-0.18.3-SNAPSHOT.jar,/home/hadoop/cloud/adam/lib/adam-core_2.10-0.18.3-SNAPSHOT.jar,/home/hadoop/cloud/adam/xubo/data/GRCH38Sub/cs-bwamem/BWAMEMSparkAll/gcdss-cli-0.0.3-SNAPSHOT.jar \
--executor-memory 4096M \
--total-executor-cores 20 BWAMEMSparkAll.jar \
hdfs://219.219.220.149:9000/xubo/callVariant/vcf/All_20160407.vcf 20

3.结果:
GRCH38chr20:

hadoop@Master:~/xubo/project/vcf/extractGRCH38chr20vcf$ ./run.sh 
start:
3360948                                                                         
load time:232552 ms
save time:243793 ms                                                             
run time:476345 ms

GRCH38chr22:

hadoop@Master:~/xubo/project/vcf/extractGRCH38chr20vcf$ ./GRCH38chr22.sh 
start:
2077087                                                                         
load time:224112 ms
save time:220048 ms                                                             
run time:444160 ms

参考

【1】https://github.com/xubo245/AdamLearning
【2】https://github.com/bigdatagenomics/adam/ 
【3】https://github.com/xubo245/SparkLearning
【4】http://spark.apache.org
【5】http://stackoverflow.com/questions/28166667/how-to-pass-d-parameter-or-environment-variable-to-spark-job  
【6】http://stackoverflow.com/questions/28840438/how-to-override-sparks-log4j-properties-per-driver

研究成果:

【1】 [BIBM] Bo Xu, Changlong Li, Hang Zhuang, Jiali Wang, Qingfeng Wang, Chao Wang, and Xuehai Zhou, "Distributed Gene Clinical Decision Support System Based on Cloud Computing", in IEEE International Conference on Bioinformatics and Biomedicine. (BIBM 2017, CCF B)
【2】 [IEEE CLOUD] Bo Xu, Changlong Li, Hang Zhuang, Jiali Wang, Qingfeng Wang, Xuehai Zhou. Efficient Distributed Smith-Waterman Algorithm Based on Apache Spark (CLOUD 2017, CCF-C).
【3】 [CCGrid] Bo Xu, Changlong Li, Hang Zhuang, Jiali Wang, Qingfeng Wang, Jinhong Zhou, Xuehai Zhou. DSA: Scalable Distributed Sequence Alignment System Using SIMD Instructions. (CCGrid 2017, CCF-C).
【4】more: https://github.com/xubo245/Publications

Help

If you have any questions or suggestions, please write it in the issue of this project or send an e-mail to me: xubo245@mail.ustc.edu.cn
Wechat: xu601450868
QQ: 601450868
使用GATK软件将vcf文件GT为1/1的纯合突变从vcf文件提取出来,可以按照以下步骤进行操作: 1. 安装并配置GATK软件:确保已经安装了GATK软件,并且已经配置了相关的环境变量和路径。 2. 读取vcf文件:使用GATKVCFReader工具读取要处理的vcf文件。 3. 过滤GT为1/1的纯合突变:使用GATK的过滤器(Filter)对vcf文件的记录进行过滤,只保留GT为1/1的纯合突变。可以使用HaplotypeCaller过滤器进行过滤。 4. 输出结果:使用GATKVCFWriter工具将过滤后的纯合突变写入新的vcf文件。 下面是使用GATK提取GT为1/1的纯合突变的步骤示例: 1. 使用命令行进入vcf文件所在的目录,并执行以下命令: ```bash gatk VCFReader -I input.vcf -O input_filtered.vcf ``` 上述命令将读取名为input.vcfvcf文件,并将其内容写入名为input_filtered.vcf的新文件。 2. 使用HaplotypeCaller过滤器对vcf文件进行过滤,只保留GT为1/1的纯合突变。执行以下命令: ```css gatk HaplotypeCaller -R reference.fasta -I input_filtered.vcf -O output_filtered.vcf --java-options "-Xmx4g" --filter "SelectNonRefSamplesWithGT==1" ``` 上述命令,reference.fasta是参考基因文件,output_filtered.vcf是过滤后的纯合突变结果文件。 3. 使用VCFWriter工具将过滤后的纯合突变写入新的vcf文件。执行以下命令: ```css gatk VCFWriter -V output_filtered.vcf -O output_cleaned.vcf ``` 上述命令将过滤后的纯合突变写入名为output_cleaned.vcf的新文件。 完成上述步骤后,output_cleaned.vcf文件就包含了GT为1/1的纯合突变了。请注意,这只是一个简单的示例,实际操作可能需要根据具体情况进行调整和优化。
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值