Нужна помощь с обнаружением пикового сигнала в Perl - PullRequest
3 голосов
/ 23 августа 2010

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

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

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.pngFrom исследование

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

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

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

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

Ответы [ 3 ]

6 голосов
/ 23 августа 2010

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

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  ==>  ..
2 голосов
/ 24 августа 2010

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

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

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

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

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

2 голосов
/ 23 августа 2010
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. Я не уверен, как вы различаете, какие из них вы хотите включить.

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...