Perl Bio::DB::Sam аварийно завершает работу во время вызова get_features_by_id()

Я использую Bio::DB::Sam в среде Centos7, используя версию 0.1.17 samtools. Я использую эту процедуру для выполнения моей установки:

wget http://sourceforge.net/projects/samtools/files/samtools/0.1.17/samtools-0.1.17.tar.bz2
tar xjf samtools-0.1.17.tar.bz2 && cd samtools-0.1.17
make CFLAGS=-fPIC
export SAMTOOLS=`pwd`
cpanm Bio::DB::Sam

который я обнаружил здесь (обратите внимание, что я изменил версию samtools)

Сбой происходит периодически, иногда на одних и тех же входных файлах. Моя общая процедура выглядит следующим образом:

  1. Используйте бабочку для создания файла .sam из файла .fastq с использованием пользовательского индекса бабочки.
  2. Используйте samtools для преобразования моего .sam в .bam, попутно сортируя и индексируя файл.
  3. Введите следующие команды Perl:

Перл:

my $sortbam = align_and_sort_and_index($reads_file);   # steps 1 and 2
my @all_gene_ids = qw(gene_id1 gene_id2 gene_id3);   # really lots more
for (my $worker=0; $worker <= $n_threads; $worker++) {
    my $pid = fork;
    die "fork error: $!" unless defined $pid;
    next if $pid;     # parent
    my @gene_ids = get_unique_subset(@all_gene_ids, $worker);
    my $sam = Bio::DB::Sam->new(-bam=>$sortbam, -fasta=>$ampl_seqfile, -autoindex=>0);
    foreach my $gene_id (@gene_ids) {
        # THIS NEXT LINE IS THE ONE THAT SEGFAULTS (SOMETIMES):
        my @alignments = $sam->get_features_by_location(-seq_id => $gene_id);
        # do something interesting with @alignments...
    }
    exit;
}

while ((my $pid=wait()) != -1) {
    print "reaped $pid\n";
}

На сегодняшний день я пробовал следующее:

  1. Увеличено количество разрешенных открытых файлов (ulimit -n)
  2. Увеличено количество разрешенных подпроцессов
  3. Увеличен лимит буферов каналов
  4. Увеличено пространство подкачки

Мы будем очень признательны за любые предложения. Благодарю вас!


person phonybone    schedule 14.06.2017    source источник
comment
Если вы запустите это без какого-либо кода fork, вы увидите ту же проблему? Если это так, это могут быть некоторые идентификаторы генов или файлы BAM, которые искажены. Я бы сначала упростил код, потому что это может быть не модуль.   -  person SES    schedule 14.06.2017
comment
Кроме того, вы упоминаете 3 разные версии samtools. Пожалуйста, уточните, что это опечатка, и один из них установлен/связан правильно.   -  person SES    schedule 14.06.2017
comment
Не имеет отношения, но взгляните на metacpan.org/pod/Parallel::ForkManager. Это сделает вашу жизнь проще.   -  person simbabque    schedule 14.06.2017
comment
Исправлены номера версий @SES samtools, спасибо.   -  person phonybone    schedule 15.06.2017
comment
@SES: я написал код, чтобы он не разветвлялся, и у меня не возникает той же проблемы. Кажется, это работает нормально. Однако код, который я изначально унаследовал, использует форк, и, как ни странно, он работает нормально при установке на машины MacOSX. Кроме того, если я вручную ограничиваю количество используемых вилок, это также работает, что изначально заставило меня заподозрить проблему с ресурсами (открытые файловые дескрипторы и т. д.).   -  person phonybone    schedule 15.06.2017
comment
@SES: рассматриваемый файл .bam создается один раз и используется в качестве входных данных для Bio::DB::Sam-›new(), поэтому я не думаю, что он искажен. Кроме того, я могу запустить версию кода с одним ответвлением для того же файла .bam, и все выглядит нормально даже после сбоя, поэтому я думаю, что с файлом все в порядке. Я еще раз проверю. Спасибо!   -  person phonybone    schedule 15.06.2017
comment
@simbabque Если я соберусь и перепишу код, я обязательно посмотрю на Parallel::ForkManager. Спасибо!   -  person phonybone    schedule 15.06.2017
comment
Что там делает оператор exit?   -  person flies    schedule 13.10.2017
comment
@flies: оператор выхода существует, потому что каждая итерация цикла порождает новый процесс (через вызов fork), затем выполняет свою работу, а затем завершает работу.   -  person phonybone    schedule 16.10.2017
comment
но это внутри цикла по идентификаторам генов, поэтому кажется, что все разветвленные процессы завершатся после просмотра первого гена? (Я почти не использовал вилку, так что IDK, если я что-то упустил.)   -  person flies    schedule 18.10.2017
comment
Я понимаю, что вы говорите; в цикле for отсутствует фигурная скобка ( '}' ). Вы правы, что он выйдет после одной итерации. Это была опечатка.   -  person phonybone    schedule 19.10.2017