Soru:
BED dosyasından çapraz hizalamalar nasıl filtrelenir?
SmallChess
2017-05-19 10:49:47 UTC
view on stackexchange narkive permalink

Bir BAM dosyam var:

  @SQ SN: chr1 LN: 248956422 @ SQ SN: chrx LN: 248956423ST-E00110: 348: HGVKKALXX: 1: 1201: 5822: 48670 323 CHR1 9999 0 67H66M16H 1000 0 GATAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAAC JJJJJJJJJJJJJJJJAJJJJJJJJJJJJFJJJJJJFJFJJJJJJFJJJJJJJJJJJA77FJFJJJ NM chrx: i: 0, MD: Z: 66: i: 66 XS: i: 65 SA: Z: chr5,18606834, - 73S76M, 34,0; RG: Z: g1  

chr1 ile hizalanmış bir okuma var ve chrx ile eşleştirilmiş.

BED dosyam var:

  chr1 0 100000 TestOnly  

BED’imin dışında kalan her şeyi filtrelemek istiyorum çapraz hizalamaları içeren bölge. Örneğimde, okumam chr1 ile hizalı olmasına rağmen, mate değil. Bunun okunmasını istemiyorum.

Okunduğumda:

samtools view -L test.bed test.bam

komut bana okumayı veriyor çünkü çapraz hizalamaları kontrol etmiyor.

Çözümüm:

samtools view -L test.bed test.bam | grep -v chrx

ancak bu çok yavaş ve sakar. Üretim hattımda şöyle bir şey yapmam gerekiyordu:

samtools view -L test.bed test.bam | grep -v chrx | grep -v ... | grep -v ... | grep -v ... | grep -v ...

S: Daha iyi bir çözüm var mı?

Bir cevap:
#1
+6
terdon
2017-05-19 22:44:29 UTC
view on stackexchange narkive permalink

SAM spesifikasyonuna göre, SAM satırının ( RNAME ) 3. alanı şudur:

RNAME: Referans dizisi NAME hizalamanın. @SQ başlık satırları mevcutsa, RNAME ("*" değilse) SQ-SN etiketlerinin birinde mevcut olmalıdır. Koordinatı olmayan, eşlenmemiş bir segment bu alanda "*" işaretine sahiptir. Bununla birlikte, eşlenmemiş bir parça, aynı zamanda, sınıflandırmadan sonra istenen bir konuma yerleştirilebilecek şekilde sıradan bir koordinata da sahip olabilir. RNAME '*' ise, POS ve CIGAR hakkında hiçbir varsayımda bulunulamaz.

Ve 7. alan (vurgu benim, "onların" yerine "eksik):

RNEXT: Şablonda okunan SONRAKİ'nin birincil hizalamasının referans sıra adı. Son okuma için, bir sonraki okuma, şablondaki ilk okumadır. @SQ başlık satırları mevcutsa, SQ-SN etiketlerinin birinde RNEXT ("*" veya "=" değilse) bulunmalıdır. Bu alan, bilgi bulunmadığında "*", ve RNEXT aynı RNAME ise "=" olarak ayarlanır . "=" Değilse ve şablondaki sonraki okumanın bir birincil eşlemesi varsa (ayrıca FLAG'deki bit 0x100'e bakın), bu alan bir sonraki okumanın birincil satırındaki RNAME ile aynıdır. RNEXT '*' ise, PNEXT ve bit 0x20 üzerinde varsayım yapılamaz

Dolayısıyla, 7. alanı = olmayan satırları kaldırmak istersiniz ve her ihtimale karşı, 7. alanı = ve 3. alan ile aynı olmayan satırlar. Bu nedenle şuna benzer bir şey kullanabilirsiniz:

  samtools view -L test.bed test.bam | awk '$ 7 == "=" || $ 3 == $ 7  

Ve tekrar bam dosyası olarak kaydetmek için:

  samtools view -L test.bed test.bam | awk '$ 7 == "=" && $ 3 == $ 7 | samtolls görünümü -b > fixed.bam  

Ayrı bir not olarak, çok nadiren birden çok grep komutunu zincirleme ihtiyaç duyulur. Bunları ayırmak için \ | (veya | 'i -E veya -P seçenekleriyle) kullanabilirsiniz. Şuna benzer bir şey:

  samtools view -L test.bed test.bam | grep -v 'chrx \ | chr2 \ | chr10 \ | chrN'  

Veya

  samtools view -L test.bed test.bam | grep -Ev 'chrx | chr2 | chr10 | chrN'  
Bunu böyle yaparsanız, 'fixed.bam` dosyasında başlık eksiktir ve bu benim deneyimime göre pek çok sorun yaratır. Başlığı her zaman tekrar eklemenizi öneririm; ya orijinal BAM'yi okurken "-h" belirterek veya ayrı olarak ekleyerek: "(samtools view -H infile.bam; samtools view…)> samtools view -b> outfile.bam`.


Bu Soru-Cevap, otomatik olarak İngilizce dilinden çevrilmiştir.Orijinal içerik, dağıtıldığı cc by-sa 3.0 lisansı için teşekkür ettiğimiz stackexchange'ta mevcuttur.
Loading...