Esistono diversi modi. Eccone due:
-
samtools idxstats ‹bamfile›
-
samtools view -H ‹bamfile› | grep '^ @ SQ'
Entrambi i comandi forniscono un elenco dei riferimenti nel file BAM. Per ottenere il loro numero, aggiungi semplicemente | wc -l
. Ma tieni presente che idxstats
restituisce anche una voce finale, *
, che non corrisponde a nessun riferimento ma a tutte le letture non mappate. Quindi sottrai questo dal numero.
Tuttavia, fai attenzione che mentre questo elenca tutti i riferimenti associati a un file BAM, non tutti i riferimenti hanno necessariamente dei record ad esso associati. Quindi potrebbero non esserci letture mappate a un dato riferimento. Questa è la terza colonna nell'output di idxstats
. Quindi, per ottenere solo quei riferimenti con almeno una lettura mappata, possiamo filtrare per quelli:
samtools idxstats ‹bamfile› | awk '$ 4! = 0' | head -n-1 | wc -l
Oh, e questi (e altri) comandi funzionano sui file BAM , non su SAM. Il comando idxstats
necessita inoltre dell'indice BAM, poiché in realtà sta leggendo le informazioni da quel file indice, non dal BAM (questo significa anche che funziona solo su file BAM ordinati). Lavorare su SAM richiederebbe il loop dell'intero file, il che è abbastanza inefficiente. In genere dovresti lavorare con BAM anziché con SAM, quindi si spera che questa non sia una restrizione.
Per fare lo stesso su un file SAM, devi iterare sul file:
awk '{print $ 3}' ‹samfile› | sort -u | wc -l
Oppure:
awk '{refs [$ 3] = 1} END {for (ref in refs) print ref}' ‹samfile ›| wc -l