dette :

  sengetøy krysser -a <myvcf>.vcf.gz -b <myinterval>.bed -wa | \ java -Xmx10g -jar snpSift.jar filter - sett <myrsid>.txt "ID i SET [0]"  

kan erstattes med en GATK SelectVariants https: // programvare .broadinstitute.org / gatk / documentation / tooldocs / current / org_broadinstitute_gatk_tools_walkers_variantutils_SelectVariants.php

  java -jar gatk.jar -T SelectVariants -R ref.fa \ -L myinterval.bed --keepIDs myrsid.txt -V myvcf.vcf.gz -o out.vcf  

og dessuten kan du også jobbe parallelt ved å dele sengefilen din (f.eks. av kromosomer), og slå sammen resultatene.

Generelt, med BEDOPS:

  $ vcf2bed < < (gunzip -c snps.vcf) | bedops -e 1 - myRegions.bed > answer.bed  

Eller:

  $ vcf2bed < < (gunzip -c snps.vcf) | bedmap --echo --echo-map-id --delim '\ t' myRegions.bed - > answer.bed  

Etc.

BEDOPS støtter inngangsstrømmer ; bruk dem der du kan.

Når det er sagt, hvis du skal spørre om VCF-dataene ofte, kan du konvertere en gang og bruke dem i nedstrømsoperasjoner.

  $ vcf2bed < < (gunzip -c snps.vcf) > snps.bed  

Deretter:

  $ bedops -e 1 snps.bed myRegions.bed > svar. seng  

Etc.

Hvis du har en beregningsklynge, kan du ytterligere parallellisere arbeidet ved å dele SNPs på kromosomer med bedextract :

  $ for chr i `bedextract --list-chr snps.bed`; gjør ekko $ {chr}; bedextract $ {chr} snps.bed > snps. $ {chr} .bed; ferdig  

Send deretter inn en batchjobb for hvert kromosom:

  $ bedmap --chrom $ {chr} --echo --echo-map-id --delim '\ t' myRegions.bed snps. $ {chr}. bed > svar. $ {chr}. bed  

Bruk bedextract og legg til - krom $ {chr} til sengekartet -anropet ditt vil gi store hastighetsforbedringer i forhold til andre verktøy, da disse alternativene fokuserer på kromosomet av interesse og ingenting annet.

Jeg fiklet med kommandoen og var i stand til å fullføre utførelsen ved å invertere rekkefølgen på rsid og sengintervallfiltrering. Kommandoen er som følger:

  gzcat <myvcf>.vcf.gz | \ java -Xmx10g -jar snpSift.jar filter - sett <myrsid>.txt "ID i SET [0]" | \ java -Xmx10g -jar snpSift.jar intervaller <mybed>.bed  

Kanskje er min rsid-delmengde mindre enn BED-intervallene mine, og dermed er prosessen raskere og mer effektiv.

Uansett takk alle

Her er en annen måte å gjøre det i 2019:

  $ bcftools view -T regions.bed -i 'ID=@<myrsid>.txt' input.vcf > output.vcf  kode> 

Ting kan økes hvis vi komprimerer og indekserer input.vcf . Da vil bcftools ha umiddelbar tilfeldig tilgang til posisjonen.

   $ bgzip -c input.vcf > input. vcf.gz $  tabix input.vcf.gz $ bcftools view -T regions.bed -i 'ID=@<myrsid>.txt' input.vcf.gz > output.vcf  

Internt fungerer bcftools med bcf . Når input og output er vcf konverterer den vcf>bcf> vcf . Når du arbeider med store filer, gir dette en stor overhead. Å konvertere vcf til bcf først kan redusere kjøretiden for videre nedstrømsanalyser dramatisk. Etter min erfaring 2-3 ganger raskere.

   $ bcftools view -Ob input.vcf > input.bcf $  bcftools index input.bcf $ bcftools view -Ob -T regions.bed -i 'ID=@<myrsid>.txt' input.bcf > output.bcf