23 Şubat 2017 Perşembe

Yeni Nesil Sekanslama Temel Veri Analizi - 2

Ufak bir uyarıyla başlayalım: YNS analizi karmaşık aşamalar içerir ancak çoğu zaman herhangi bir aşamanın tam olarak nasıl kullanılacağını bilmenizi gerektirmez. Bu nedenle bu ve sonraki yazıları okurken anlamadığınız bir kısım olursa tekrar okuyun, eğer yine de karışık geliyorsa sadece ilgili komutları çalıştırıp devam edin. Tüm süreci genel olarak anladığınızı düşündüğünüzde ilgili kısma tekrar göz atıp en azından ne amaçla kullanıldığını genel hatlarıyla anlamaya çalışın.

Yeni Nesil Sekanslama deneylerindeki en önemli aşamalardan birisi, kütüphane hazırlama aşamasıdır. Yaygın olarak kullanılan Illumina teknolojisinde sekanslanması amaçlanan DNA bölgesi -kullanılacak sarf malzemesinin özelliklerine göre- kısa DNA parçalarına bölünür ve iki ucuna adaptörler takılır. Adaptörler, sekanslanacak DNA parçasının sekanslama yüzeyindeki oligonükleotitlere hibridize olması için kullanılan ve dizilimi özel olarak belirlenen oligonükleotitlerdir. Kütüphane hazırlama aşaması tamamlandığında her bir DNA parçası aşağıdaki gibi bir yapıya kavuşur:




5' ucuna eklenen adaptörün bittiği yerden itibaren okuma başlar ve kullanılan sarf malzemesinin kapasitesi kadar baz okunur. Bu durumda ya okunabilecek baz sayısı sekanslanacak DNA parçasından küçüktür, ya DNA parçasının uzunluğuna eşittir, ya da DNA parçasının uzunluğundan büyüktür. Bahsettiğim son durumda okumanın içerisine 3' ucundaki adaptörün bir kısmı da dahil olur. Örneğin, indirdiğiniz FASTQ dosyası her biri 50'şer bazlık okumalar içeriyor. 30 baz uzunluğunda bir DNA parçasınının sekanslandığı okumalarda kalan 20 baz 3' adaptörden okunuyor demektir. Bu gibi durumlarda yapılacak temizleme işlemi 3' adaptörlerinin çıkarılması ile başlar. Bu amaçla kullanılabilecek birçok araç mevcut ancak cutadapt isimli yazılımı kullanmayı tercih ediyorum. Programı kurduktan sonra ihtiyacımız olan tek şey 3' adaptörün içeriği. Veriyi aldığım makalenin ekinde adaptör içeriği belirtilmiş: CTGTAGGCACCATCAAT . cutadapt'ın kullanışlı özelliklerinden birisi de, doğrudan sıkıştırılmış  formattaki (.gz uzantılı) dosyalar üzerinden işlem yapabilmesi. Bunun için çalıştıracağımız komut şu şekilde olacak:

cutadapt -a CTGTAGGCACCATCAAT -e 0.1 -O 5 -m 15 -o  SRR1513652_trimmed.fastq SRR1513652.fastq.gz

Hangi sekansın adaptör olduğunu -a parametresi ile belirtiyorum. -e parametresi hata oranını belirtiyor: bu değeri 0.1 olarak belirtmekle her on bazda bir bazlık eşleşme hatasını kabul edebileceğimizi  ifade ediyoruz. -O, adaptörün en az kaç bazlık kısmının eşleşmesi gerektiğini belirlemek için kullanıyoruz, -m ise bir okumadan adaptör çıkarıldığında kalan kısmın en az kaç baz uzunluğunda olması gerektiğini, bu değeri geçmeyen okumaların filtrelenmesini belirtmek için kullanılıyor. Yani adaptör çıkarıldığında kalan okuma uzunluğu en az 15 baz olmalı diyoruz. Son olarak kullandığımız parametre olan -o ile temizlenmiş okumaların hangi dosyaya yazılacağını belirtiyoruz (eğer aynı klasörde aynı isme sahip başka bir dosya varsa o dosyanın üzerine yazılır). Parametreleri yazmayı tamamladıktan sonra artık sıra hangi dosyanın temizlenmesi gerektiğini yazmaya geldi, bunun için de indirmiş olduğumuz SRR1513652.fastq.gz isimli dosyayı kullanıyoruz. Terminal'de çalıştığımız klasörün içerisinde bu dosyanın yer aldığını kontrol etmeyi unutmayın, geçen yazımızda bu dosyayı Masaüstü (Desktop) klasörüne kopyalamıştık. Programı çalıştırdıktan sonra ekranda birçok şey yazmakla birlikte şöyle bir ifade yer alıyor:

=== Summary ===

Total reads processed:               2,177,511
Reads with adapters:                 2,173,988 (99.8%)
Reads that were too short:               3,686 (0.2%)
Reads written (passing filters):     2,173,825 (99.8%)

Toplam 2.177.511 okumanın işlendiği belirtiliyor, hatırlarsanız geçen yazımızda dosyanın içerisindeki okuma sayısını hesapladığımızda da aynı sayıyı elde etmiştik. Okumaların %99.8'inde adaptöre rastlandığı belirtiliyor, bunun temel sebebi de sekansladığımız DNA parçalarının en fazla 40 baz civarında olmasını beklememiz, buna ilişkin detayları makalede bulabilirsiniz. Temizlenen okumaların %0.2'si 15 bazdan daha kısa, bu alt sınırı -m parametresi ile belirtmiştik. Filtreleri geçen ve temizlenen okuma sayısını da en alt satırda görüyoruz. Ekrandaki bilgilerin detayına daha fazla girmeden oluşturulan dosyanın içeriğine bir göz atalım:

head SRR1513652_trimmed.fastq

Bu komutu çalıştırdıktan sonra beklediğimiz çıktı şöyle olmalı:

@SRR1513652.1 HWI-ST570:67:D10FDACXX:4:1101:1428:2373 length=50
ACTGGTCTTCATTAATCAAGAACGAAAGTCA
+
@CCFFFFFHHDHHJIIGHHGI:FFCGIIHHG
@SRR1513652.2 HWI-ST570:67:D10FDACXX:4:1101:2008:2430 length=50
ACTGCACCGTGACTCTCCTCGAAGACAGTGTCAAGCGG
+
CCCFFFFFHHHHHJIJIIJJJ1CGHIJIIHHIJJGIIJ
@SRR1513652.3 HWI-ST570:67:D10FDACXX:4:1101:2343:2302 length=50
ACTGTGTTCAGTTTATCGTACAAAT

Her bir okumanın kısaldığını görüyorsunuz değil mi? Şimdi sırada sadece DNA dizilimlerini almak var, bunun için yeni bir komut kullanacağız: paste. Bu komut, alt alta yer alan satırların belirli bir motif takip ederek yan yana yazılması için kullanılıyor. Bu komuta neden ihtiyacımız olduğunu bir sonraki komutu kullandığımızda anlayacağız, şimdilik ne işe yaradığını görelim. FASTQ dosyasındaki her bir okuma 4 satırlık parçalar halinde yer alıyor, her bir okumayı bir satırda görüntülemek için paste komutunu şu şekilde kullanmalıyız:

paste - - - - < SRR1513652_trimmed.fastq > SRR1513652.paste 

Yukarıda yer alan ve aralarında boşluk bulunan dört adet - karakteri, ardarda gelen her dört satırın yan yana yazılması gerektiğini söylüyor. < karakterinin özel bir anlamı var: paste komutuna girdi olarak kullanılacak dosyanın ne olduğunu belirtiyor. < karakterini, ilgili komutu dosyayla besleyen veya dosyayı ilgili komuta yönlendiren bir huni gibi düşünebilirsiniz. Benzer şekilde > karakteri de elde edilen çıktıların bir dosyaya yazılmasını/yönlendirilmesini sağlıyor, bu durumda elde ettiğimiz çıktıları SRR1513652.paste dosyasına yazdırıyoruz. Yukarıdaki komutu çalıştırıp elde edilen dosyaya da 

head SRR1513652.paste

komutunu kullanarak baktığınızda her bir okumanın yan yana ve tab ile ayrılmış olarak görüntülendiğini göreceksiniz. Şimdi yapmamız gereken şey, ikinci sütunları alıp bunları bir dosyaya yazmak. Bunun için de yeni bir komutumuz var :) Sürekli yeni komutlara maruz kaldığınız için yorucu bir süreç olduğunun farkındayım, sabredin, az kaldı :) Burada kullanacağımız komut, her derde deva awk komutu. awk, metin işleme, değiştirme ve dönüştürme için hem çok kullanışlı, hem de çok hızlı bir UNIX programı. Burada dikkat etmenizi istediğim şey, awk'ın her bir satırı ayrı ayrı okuyarak üzerinde istediğimiz işlemi yaptığıdır. Tab ile ayrılmış dosyadaki her bir satırdaki sadece ikinci sütunu almak ve bunu sonrasında bir dosyaya yazmak için şu komutu kullanıyoruz:

awk -F\t '{print $2}' SRR1513652.paste > SRR1513652.seq

Ve işte bu kadar! Artık elimizde her bir okumaya ait DNA dizilimleri var! Yukarıdaki satırda ne yaptığımızı özetleyerek başlayayım. -F parametresi, okunan dosyadaki her bir satırın sütunlara dönüştürülürken hangi karakter kullanılacağını belirtmek için kullanılır. Buradaki sütunları, bir Excel dosyasındaki sütunlar olarak düşünebilirsiniz. Yukarıda kullandığımız paste komutu her bir dörtlü satırı tab karakteriyle birleştirdiği için, bu dörtlü satırları tekrar sütunlara bölmek için tab karakterini kullanıyoruz. Tab karakteri bilgisayar dünyasında \t ile ifade edilir. Burada dikkatinizi çekmek istediğim bir şey var: Şimdiye kadar kullandığımız komut veya programlarda parametleri belirtirken aralarda mutlaka birer boşluk bıraktık ancak -F parametresi bir istisna oluşturuyor ve ayıraç olarak kullandığımız \t karakteri -F'ye birleşik olarak kullanılıyor. awk kullanarak ne yapmak istediğimizi '{ ve }' arasına koyarak belirtiyoruz. awk'ı bundan sonra sıklıkla kullanacağımız için bu tür detaylara dikkat etmenizi öneririm. Herhangi bir ifadeyi ekrana yazdırmak için awk'a özel print fonksiyonunu kullanıyoruz ve arada bir boşluk bırakarak ekrana neyi yazdırmak istediğimizi belirtiyoruz. Birinci sütun $1 isimli değişkenin içerisinde yer alıyor, ikincisi $2'de, üçüncüsü $3'te, ve bu böyle gidiyor. Okumamız ikinci sütunda yer aldığı için biz ekrana sadece ilgili sütünü yazdırıyoruz. Hangi dosya üzerinde işlem yapacağımızı, veya hangi dosyayı satır satır okuyacağımızı belirtmek için bir boşluğun ardından dosyamızın adını yazıyoruz. Normal şartlarda ekrana yazdıracağımız bilgileri de hunimiz ( ) yardımıyla SRR1513652.seq isimli dosyaya yönlendiriyoruz. Bu aşamada SRR1513652.seq isimli bir dosya oluşturuluyor ve ekrana yazdırılacak çıktılar doğrudan bu dosyaya yazdırılıyor. Gelelim bu dosyanın içeriğine:

head SRR1513652.seq

komutunu çalıştırdığımızda her şeyin yolunda olduğunu görebiliriz:

ACTGGTCTTCATTAATCAAGAACGAAAGTCA
ACTGCACCGTGACTCTCCTCGAAGACAGTGTCAAGCGG
ACTGTGTTCAGTTTATCGTACAAAT
ACTGTTCTTAGTTGGTGGAGTGATTTGTCTGGTTT
ACTGTATGAGATAGTGACTAAAAATATAAA
ACTGTTTTATTTGGAGCATGATCAA
ACTGATTCCACTTGGGGTAAATTC
ACTGGACCATTTCGTCCAGAGTTTTTGGGT
ACTGTCAACTGGAAATATTTATCAA
ACTGCACCTTATCCGGGATCCTCATATGG

Şimdiye kadarki süreçte birçok dosya oluşturduk ve ortalık biraz karıştı. Elimizdeki sekanslarla oynamaya başlamadan önce çalıştırdığımız tüm komutları derli toplu bir hale getirip daha sonra rahatlıkla kullanabileceğimiz tek satırlık bir komut dizisine dönüştürsek fena olmayacak, hem böylece bu komut dizisini birden çok dosya ile çalışırken de rahatlıkla kullanabiliriz. 

Bir sonraki yazıda bu işlemin nasıl yapılacağını anlatacağım ve elimizdeki okumalara ilişkin temel istatistikleri nasıl elde edebileceğimizi anlatmaya başlayacağım.

Bu bölümde kullandığımız program: cutadapt
Bu bölümde kullandığımız UNIX komutları: head, paste, awk, > ve < .
Bu bölümde kullandığımız dosyalar: SRR1513652.fastq.gz, SRR1513652_trimmed.fastq, SRR1513652.paste ve SRR1513652.seq .