Извлеките определенные столбцы и строки из файла, используя Linux или Python - PullRequest
Новая
       7

Извлеките определенные столбцы и строки из файла, используя Linux или Python

0 голосов
/ 28 июня 2018

Я столкнулся с проблемой при работе с файлом 12 Гб. Я новичок, когда дело доходит до Linux. Я надеюсь, что кто-то здесь может помочь мне, любые предложения приветствуются.

У меня есть файл с именем phase_3.vcf , который выглядит так: 

##INFO=<ID=EAS_AF,Number=A,Type=Float,Description="Allele frequency in the EAS populations">                            
##INFO=<ID=EUR_AF,Number=A,Type=Float,Description="Allele frequency in the EUR populations">                            
##INFO=<ID=AFR_AF,Number=A,Type=Float,Description="Allele frequency in the AFR populations">                            
##INFO=<ID=AMR_AF,Number=A,Type=Float,Description="Allele frequency in the AMR populations">                            
##INFO=<ID=SAS_AF,Number=A,Type=Float,Description="Allele frequency in the SAS populations">                
#CHROM  POS ID  REF ALT QUAL    FILTER  INFO
1   10177   rs367896724 A   AC  .   .    dbSNP_150;E_Freq;E_1000G;EAS_AF=0.3363;SAS_AF=0.4949;AFR_AF=0.4909
1   10235   rs540431307 T   TA  .   .   dbSNP_150;E_Freq;E_1000G;EAS_AF=0.0000;AMR_AF=0.0014;
1   10352   rs555500075 T   TA  .   .   dbSNP_150;EAS_AF=0.4306;EUR_AF=0.4264;SAS_AF=0.4192;
1   10505   rs548419688 A   T   .   .   TSA=SNV;E_Freq;EAS_AF=0.0000;AMR_AF=0.0000;AFR_AF=0.0008
1   10506   rs568405545 C   G   .   .   dbSNP_150;TSA=SNV;MA=G;MAF=0.000199681;EAS_AF=0.0000;

И я хочу сохранить первые 5 столбцов, строки "EAS_AF =" и цифры, следующие за ними.

Для простоты ожидаемая форма результата с именем result.txt должна быть такой: 

#CHROM  POS ID  REF ALT INFO
1   10177   rs367896724 A   AC  EAS_AF=0.3363
1   10235   rs540431307 T   TA  EAS_AF=0.0000
1   10352   rs555500075 T   TA  EAS_AF=0.4306
1   10505   rs548419688 A   T   EAS_AF=0.0000
1   10511   rs534229142 G   A   EAS_AF=0.0000
1   10539   rs537182016 C   A   EAS_AF=0.0000

Ответы [ 2 ]

0 голосов
/ 28 июня 2018

Ваш файл в очень распространенном формате биоинформатики (vcf). Итак, в биоинформатике есть инструменты, специально предназначенные для решения вашей задачи: Например, bcftools имеет возможность удалить все элементы INFO, кроме одного.

Недостатком этого инструмента является то, что он строго требует vcf-формат. Так что это приведет к ошибкам на вашем примере. Но я думаю, что вы сократили заголовок для этого поста, и он должен быть в порядке в исходном файле. Чтобы сделать его полезным для меня, мне пришлось настроить заголовок, как описано в определении формата , добавив формат файла и строку INFO в заголовке для каждой другой информации, которую вы аннотировали в своих вариантах: 

##fileformat=VCFv4.1
##FILTER=<ID=PASS,Description="All filters passed">
##INFO=<ID=dbSNP_150,Number=0,Type=Flag,Description="has dbSNP 150 entry">
##INFO=<ID=E_Freq,Number=0,Type=Flag,Description="has E_Freq entry">
##INFO=<ID=E_1000G,Number=0,Type=Flag,Description="has E_1000G entry">
##INFO=<ID=EAS_AF,Number=A,Type=Float,Description="Allele frequency in the EAS populations">
##INFO=<ID=EUR_AF,Number=A,Type=Float,Description="Allele frequency in the EUR populations">
##INFO=<ID=AFR_AF,Number=A,Type=Float,Description="Allele frequency in the AFR populations">
##INFO=<ID=AMR_AF,Number=A,Type=Float,Description="Allele frequency in the AMR populations">
##INFO=<ID=SAS_AF,Number=A,Type=Float,Description="Allele frequency in the SAS populations">
##INFO=<ID=TSA,Number=1,Type=String,Description="No idea">
##INFO=<ID=MA,Number=1,Type=String,Description="Minor Allele">
##INFO=<ID=MAF,Number=1,Type=Float,Description="Minor Allele Frequency">
##contig=<ID=1>
##bcftools_annotateVersion=1.4-23-ga468a83+htslib-1.4-34-g8e1be4a
##bcftools_annotateCommand=annotate -x QUAL test2.vcf.gz; Date=Thu Jun 28 17:58:17 2018
#CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO
1       10177   rs367896724     A       AC      .       .       dbSNP_150;E_Freq;E_1000G;EAS_AF=0.3363;SAS_AF=0.4949;AFR_AF=0.4909
1       10235   rs540431307     T       TA      .       .       dbSNP_150;E_Freq;E_1000G;EAS_AF=0;AMR_AF=0.0014
1       10352   rs555500075     T       TA      .       .       dbSNP_150;EAS_AF=0.4306;EUR_AF=0.4264;SAS_AF=0.4192
1       10505   rs548419688     A       T       .       .       TSA=SNV;E_Freq;EAS_AF=0;AMR_AF=0;AFR_AF=0.0008
1       10506   rs568405545     C       G       .       .       dbSNP_150;TSA=SNV;MA=G;MAF=0.000199681;EAS_AF=0

В зависимости от вашего заголовка вам может потребоваться сделать два дополнительных шага, чтобы использовать этот инструмент. Во-первых, это сжать ваш ввод с помощью bgzip: 

bgzip phase_3.vcf

Второе - создание tabix-индекса для быстрого доступа к вашему архивному файлу (в качестве выходного файла создается дополнительный файл phase_3.vcf.gz.tbi): 

tabix phase_3.vcf.gz

Фактический вызов bcftools после того, как вы ввели в нужном формате, просто: 

bcftools annotate -x ^INFO/EAS_AF phase_3.vcf.gz >phase_3_output.vcf

С помощью этих шагов я получаю что-то очень близкое к желаемому результату: 

##fileformat=VCFv4.1
##FILTER=<ID=PASS,Description="All filters passed">
##INFO=<ID=EAS_AF,Number=A,Type=Float,Description="Allele frequency in the EAS populations">
##contig=<ID=1>
##bcftools_annotateVersion=1.4-23-ga468a83+htslib-1.4-34-g8e1be4a
##bcftools_annotateCommand=annotate -x ^INFO/EAS_AF test2.vcf.gz; Date=Thu Jun 28 17:59:04 2018
#CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO
1       10177   rs367896724     A       AC      .       .       EAS_AF=0.3363
1       10235   rs540431307     T       TA      .       .       EAS_AF=0
1       10352   rs555500075     T       TA      .       .       EAS_AF=0.4306
1       10505   rs548419688     A       T       .       .       EAS_AF=0
1       10506   rs568405545     C       G       .       .       EAS_AF=0

Удаление первых нескольких строк с заголовком, удаление столбцов QUAL и FILTER с обрезкой и перезапись от 0 до 0,0000 с помощью sed завершает вашу задачу: 

bcftools annotate -x ^INFO/EAS_AF phase_3.vcf.gz | tail -n +7 | cut -f1-5,8 |sed 's/=0$/=0.0000/g' >phase_3_output_finished.vcf

И вы получите желаемый результат: 

#CHROM  POS ID  REF ALT INFO
1   10177   rs367896724 A   AC  EAS_AF=0.3363
1   10235   rs540431307 T   TA  EAS_AF=0.0000
1   10352   rs555500075 T   TA  EAS_AF=0.4306
1   10505   rs548419688 A   T   EAS_AF=0.0000
1   10511   rs534229142 G   A   EAS_AF=0.0000
1   10539   rs537182016 C   A   EAS_AF=0.0000

В зависимости от вашей конечной цели, вы можете даже найти способы ее достижения без этого шага, о котором говорилось здесь, узнав об опциях таких инструментов, как bcftools .

Поскольку вы работаете с данными биоинформатики, вы можете найти дополнительную помощь в соответствующих сообществах, таких как biostars или bioinformatics.stackexchange

0 голосов
/ 28 июня 2018

Избавьтесь от первых строк, если не нужно, а затем используйте pandas для импорта csv и присвойте имена столбцам dataframe 

    import pandas as pd
    import numpy as np
    df = pd.read_csv('/example.csv',names=['Column1','Column2'....])

edit: это скрипт на python, вам нужно установить python в вашей системе linux или использовать тот, который встроен в версию.

...