Bio :: DB :: Sam - Получить количество отображений для всех операций чтения в файле bam.

Question

Я хочу вычислить выражения стенограммы и, следовательно, мне нужно получить количество отображений для всех операций чтения в файле bam.Моя текущая процедура состоит в том, чтобы пройти полные расшифровки и получить чтения, которые отображены на нем, используя Bio :: DB :: Sam.Результаты сохраняются в хэше с read_name в качестве ключа (10 букв) и number_of_mappings в качестве значения (целое число).

Вот код, который я использую:

use Bio::DB:Sam;
use strict;

my %global_read_occurrences;


sub getGlobalReadOccurrences {

 my ($ids, $bam_file) = @_;

 $sam = Bio::DB::Sam -> new (-bam => $bam_file);

 foreach my $id (@{$ids}){
   my $alignments = $sam -> get_features_by_location(-seq_id => $transcript_id, -iterator => 1);


  while (my $alignment = $alignments -> next_seq){

   my $read_name = $alignment -> query -> name;

   if (exists($global_read_occurrences{$read_name})){
    $global_read_occurrences{$read_name}++;
   }
   else {
    $global_read_occurrences{$read_name} = 1;
   }
  }
 }
}

Мои вопросы:Есть ли какая-либо другая возможность, где я могу получить количество глобальных сопоставлений на чтение непосредственно и где мне не нужно просматривать все стенограммы?Я не смог найти никаких подпрограмм в Bio :: DB :: Sam, таких как $ sam -> getNumberOfMappings ($ read_name);

Я использую файлы bam с более чем 50 миллионами сопоставленных операций чтения, поэтому хеш идетнужны огромные ресурсы памяти (иногда около 40 ГБ) Это действительно возможно или это происходит откуда-то еще?И есть ли другая возможность хранить данные с меньшим количеством памяти?

Большое спасибо!

Answer 1

Файлы BAM обычно сортируются по хромосомному расположению, а не по имени чтения, поэтому сопоставления чтения могут быть расположены в любом месте файла. Самое простое, что вам нужно сделать - это зайти в файл SAM и запустить простую команду оболочки:

 cut -f1,1 myfile.sam | sort | uniq -c

Это создаст файл, подобный этому:

  2 HWI-EAS299_4_30M2BAAXX:2:99:965:826
  2 HWI-EAS299_4_30M2BAAXX:2:99:966:1932
  2 HWI-EAS299_4_30M2BAAXX:2:99:971:146
  2 HWI-EAS299_4_30M2BAAXX:2:9:997:1263
  2 HWI-EAS299_4_30M2BAAXX:2:99:972:281
  2 HWI-EAS299_4_30M2BAAXX:2:99:973:1904
  1 HWI-EAS299_4_30M2BAAXX:2:99:976:186
  2 HWI-EAS299_4_30M2BAAXX:2:99:986:687
  6 HWI-EAS299_4_30M2BAAXX:2:99:987:165
  2 HWI-EAS299_4_30M2BAAXX:2:99:99:1582
  2 HWI-EAS299_4_30M2BAAXX:2:99:99:160
  2 HWI-EAS299_4_30M2BAAXX:2:99:998:1139

Первый столбец - это количество отображений. Второе - это прочитанное имя.