Ha a klikkelős, bugos, zárt forráskódú mikroRNS elemző programokból elegünk lesz, érdemes megnézni, milyen alternatívákkal dolgozhatunk. Az egyik ezek közül a mirDeep-P. A telepítése nem igényel különleges előkészítéseket. A Bowtie legtöbb bioinformatikus kelléktárában már megtalálható, a Vienna package pedig egyszerű C programok halmaza.
Az első komolyabb lépés a szekvenciák átalakítása egy redundanciától mentes formára, ami így néz ki:
>azonosito_xszam
szekvencia
Mivel a szekvenciák elég rövidek, célszerű azt használni azonosítónak is, mert az eredemény fájlokból könnyen meglesz a szekvencia is. Ezzel az egyszerűsítéssel a szekvenciák átalakítására elég egy parancssori szkript:
awk '{if( (NR-1) % 4 == 1)print}' input.fastq | sort | uniq -c | awk '{print ">"$2"_x"$1"\n"$2}' >input.fasta
Mit is csinál ez? Az első awk kiszedi minden negyedik sort a második sortól kezdődően. Sorba rendezzük, megszámoljuk mennyi van, majd a második awk a fenti formátumra hozza ezt.
A vizsgált élőlény genomját beindexeljük Bowtie-al. Ez elég egyszerű lépés, de akadhatnak problémák. Például a búza genomja olyan nagy, hogy az index nem fér bele a 32 bites számábrázolásba. Az egyik lehetőség a bowtie2 használata, ami támogatja a 64 bites indexeket, vagy az A, B, D genomonkénti indexelés.
Szükség lesz még a potenciális r/t/s/snoRNS-ek helyére GFF formátumban. Jól annotált élőlényeknél ezt könnyű letölteni, de előfordulhat, hogy nekünk kell elkészíteni. Az RFam adatbázisból letölthetőek az ismert RNS-ek, egy Blast segítségével magunk is készíthetünk közelítő megoldást.
Miután minden feltétel teljesült, ideje futtatni a szkripteket.
bowtie -p 6 -a -v 0 bindex -f input.fasta >align_to_genome.tmp
Sok magyarázatot nem kell fűzni. A readeket illesztjük a genomra.
convert_bowtie_to_blast.pl align_to_genome.tmp input.fasta reference.fasta >align
A Bowtie kimenetet táblázattá alakítjuk.
filter_alignments.pl align -c 15 >depth.filtered
Szűrűnk lefedettség alapján.
overlap.pl depth.filtered data/ncRNA.gff -b >overlap_ncrna
Megnézzük mely illesztések fednek át más RNS típusokkal.
alignedselected.pl depth.filtered -g overlap_ncrna >noncrna
Kiszűrjük, a szükségtelen illesztéseket.
filter_alignments.pl noncrna -b input.fasta >candidate.mirna.fasta
A filter_alignments másodszori futtatásával a readeket szűrjük. Ez a kimenet majd egy későbbi lépésnél lesz fontos.
excise_candidate.pl reference.fasta noncrna 250 >precursors.fa
Meghatároztuk a prekurzorokat.
cat precursors.fa | RNAfold --noPS >hairpins.txt
Az RNS prekurzor térbeli struktúráját határozzuk meg. Ez a lépés a leglassabb és itt érdemes párhuzamosítani. Először is a prekurzorok lehetnek redundánsak, különösen egy repetitív elemekkel gazdagított növény esetén. Ha a prekurzos szekvencia nukleotidra pontosan megegyezik, felesleges kétszer lefuttatni ugyan azt az algoritmust. A következő szkript eltünteti a redundanciát:
grep -v "^>" precursors.fa | sort -u | awk '{print ">"NR"\n"$1}' >precursors.nored.fa
Ha van lehetőség, hogy farmon vagy több gépen elosztva futtassunk, akkor célszerű a precursors.nored.fa fájlt több darabra szétszedni és úgy futtatni. A térszerkezetek megfeleltetése az eredeti kandidáns mikroRNS prekurzoroknak egy kicsit nehezebb, több kódolást igényel.
bowtie-build -f precursors.fa bindex/precursors
Készítünk egy Bowtie indexet a prekurzorokból.
bowtie -a -v 0 bindex/precursors -f candidate.mirna.fasta >align_to_precursors.tmp
Emlékszünk még a kandidáns miRNS szekvenciákra? Itt az idő, hogy felhasználjuk őket és a prekurzorokra illesszük őket.
convert_bowtie_to_blast.pl align_to_precursors.tmp candidate.mirna.fasta precursors.fa >align2
Régi ismerősünk, a Bowtie kimenet konvertáló ismét akcióban.
sort +3 -25 align2 >signatures
Ez kell, azt mondja a dokumentáció. Gondolom egyszerűbb volt külső programot használni a rendezésre. Rendezett adatokra pedig hatékonyabb kódot lehet írni.
miRDP.pl signatures hairpins.txt >prediction
Végül az egész folyamat megkoronázásaként jöjjön a végső program. Növények esetén van még egy opcionális lépés, ami egy megmondó ember kritériumai alapján meghatározza, mi lehet "igazi" miRNS. Tényleg nem akarom kissebbíteni Meyers érdemeit, de bármelyik szaktekintély állíthat fel kritérium rendszert, amit majd egy következő szaktekintély félresöpörhet. A reference.len a referencia szekvencia azonosítói és a szekvenciák hossza.
rm_redundant_meet_plant.pl reference.len precursors.fa prediction nr.prediction filter.prediction
Az eredményeket végül felhasználhatjuk, hogy megírjuk a következő nagy cikkünket.