Нужна помощь с обнаружением пикового сигнала в Perl


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

Пример некоторых значений

5.7
5.3
8.2
16.5
34.2
58.8
**75.4**
75
65.9
62.6
58.6
66.4
71.4
53.5
40.5
26.8
14.2
8.6
5.9
7.7
14.9
30.5
49.9
69.1
**75.3**
69.8
58.8
57.2
56.3
67.1
69
45.1
27.6
13.4
8
5
Эти значения показывают два пика на 75.4 и 75.3, вы можете видеть, что значения увеличиваются, а затем уменьшаются. Изменение не всегда одно и то же.

График интенсивности значения

http://lh4.ggpht.com/_aEDyS6ECO8s/THKTLgDPhaI/AAAAAAAAAio/HQW7Ut-HBhA/s400/peaks.png из исследования

Одна из вещей, которые я думаю сделать, это хранить каждую из групп, т. е. горы в хэше, а затем искать наибольшее значение в группе. Один из вопросов, которые я вижу, заключается в том, как определить границы каждой из групп.

Вот ссылка на код, который у меня есть так далеко: http://paste-it.net/public/y485822/

Вот ссылка на полный набор данных: http://paste-it.net/public/ub121b4/

Я пишу свой код на Perl. Любая помощь будет очень признательна. Спасибо

3 3

3 ответа:

Вам нужно решить, насколько локальными вы хотите, чтобы пики были. Подход здесь находит пики и впадины в широких областях данных.

use strict;
use warnings;

my @data = (
    5.7, 5.3, 8.2, 16.5, 34.2, 58.8, 75.4, 75, 65.9, 62.6,
    58.6, 66.4, 71.4, 53.5, 40.5, 26.8, 14.2, 8.6, 5.9, 7.7,
    14.9, 30.5, 49.9, 69.1, 75.3, 69.8, 58.8, 57.2, 56.3, 67.1,
    69, 45.1, 27.6, 13.4, 8, 5,
);

# Determine mean. Or use Statistics::Descriptive.
my $sum;
$sum += $_ for @data;
my $mean = $sum / @data;

# Make a pass over the data to find contiguous runs of values
# that are either less than or greater than the mean. Also
# keep track of the mins and maxes within those groups.
my $group = -1;
my $gt_mean_prev = '';
my @mins_maxs;
my $i = -1;

for my $d (@data){
    $i ++;
    my $gt_mean = $d > $mean ? 1 : 0;

    unless ($gt_mean eq $gt_mean_prev){
        $gt_mean_prev = $gt_mean;
        $group ++;
        $mins_maxs[$group] = $d;
    }

    if ($gt_mean){
        $mins_maxs[$group] = $d if $d > $mins_maxs[$group];
    }
    else {
        $mins_maxs[$group] = $d if $d < $mins_maxs[$group];
    }

    $d = {
        i       => $i,
        val     => $d,
        group   => $group,
        gt_mean => $gt_mean,
    };
}

# A fun picture.
for my $d (@data){
    printf
        "%6.1f  %2d  %1s  %1d  %3s  %s\n",
        $d->{val},
        $d->{i},
        $d->{gt_mean} ? '+' : '-',
        $d->{group},
        $d->{val} == $mins_maxs[$d->{group}] ? '==>' : '',
        '.' x ($d->{val} / 2),
    ;

}

Вывод:

   5.7   0  -  0       ..
   5.3   1  -  0  ==>  ..
   8.2   2  -  0       ....
  16.5   3  -  0       ........
  34.2   4  -  0       .................
  58.8   5  +  1       .............................
  75.4   6  +  1  ==>  .....................................
  75.0   7  +  1       .....................................
  65.9   8  +  1       ................................
  62.6   9  +  1       ...............................
  58.6  10  +  1       .............................
  66.4  11  +  1       .................................
  71.4  12  +  1       ...................................
  53.5  13  +  1       ..........................
  40.5  14  -  2       ....................
  26.8  15  -  2       .............
  14.2  16  -  2       .......
   8.6  17  -  2       ....
   5.9  18  -  2  ==>  ..
   7.7  19  -  2       ...
  14.9  20  -  2       .......
  30.5  21  -  2       ...............
  49.9  22  +  3       ........................
  69.1  23  +  3       ..................................
  75.3  24  +  3  ==>  .....................................
  69.8  25  +  3       ..................................
  58.8  26  +  3       .............................
  57.2  27  +  3       ............................
  56.3  28  +  3       ............................
  67.1  29  +  3       .................................
  69.0  30  +  3       ..................................
  45.1  31  +  3       ......................
  27.6  32  -  4       .............
  13.4  33  -  4       ......
   8.0  34  -  4       ....
   5.0  35  -  4  ==>  ..
my @data = ...;

# filter out sequential duplicate values
my @orig_index = 0;
my @deduped = $data[0];
for my $index ( 1..$#data ) {
    if ( $data[$index] != $data[$index-1] ) {
        push @deduped, $data[$index];
        push @orig_index, $index;
    }
}

# add a sentinel (works for both ends)
push @deduped, -9**9**9;

my @local_maxima_indexes;
for my $index ( 0..$#deduped-1 ) {
    if ( $deduped[$index] > $deduped[$index-1] && $deduped[$index] > $deduped[$index+1] ) {
        push @local_maxima_indexes, $orig_index[$index];
    }
}
Обратите внимание, что при этом первое значение считается локальным максимумом, а также значения 71.4 и 69. Я не уверен, как вы различаете, какие из них вы хотите включить.

Есть ли у вас контрольный набор данных? Если это так, я бы рекомендовал нормализовать ваши данные, используя, скажем, простое логарифмическое соотношение между интенсивностями дрожжей и контрольными изображениями.

Затем вы можете использовать порт perl ChiPOTle для захвата значимых пиков, что звучит намного надежнее, чем поиск локальных/глобальных максимумов и т. д.

ChiPOTle-это алгоритм поиска пиков, используемый для анализа данных микрочип-микросхемы, но я успешно использовал его во многих других приложениях (например, ChIP-seq, который по общему признанию, ближе к своему первоначальному назначению, чем в вашем случае).

Полученные отрицательные значения логарифма (дрожжи / контроль) будут использованы для построения гауссовой фоновой модели оценки значимости. Затем алгоритм использует частоту ложных обнаружений для многократной коррекции тестирования.

Вот оригинал статьи .