Сглаживание с низкостью - PullRequest
       12

Сглаживание с низкостью

1 голос
/ 13 октября 2011

Я пытаюсь преобразовать скрипт Matlab в R, и у меня возникли проблемы со сглаживанием.

Код Matlab, который я хочу преобразовать, следующий:

for i = 1:size(spike_sum,2)
        smooth_sum(1:Ne,i)=smooth(double(spike_sum(1:Ne,i)),spanNe,'lowess');
end

for i = 1:Ne
        smoother_sum(i,:)=smooth(double(smooth_sum(i,:)),spanT,'lowess');
end

где spike_sum - это матрица с Ne x 4000. Я хочу сначала сгладить в Dim 1 с span spanNe, и сделать это для всех 4000 срезов. Затем я хочу сгладить в Dim 2 с span spanT, и сделать это для всех Ne срезов.

Я смотрел на функцию lowess в R, но кажется, что она принимает в двух измерениях lowess (x, y, span, iter, delta). Итак, чтобы получить результаты кода выше в R, я бы просто взял кусок матрицы для y и повторил бы постоянное значение для x?

Ответы [ 2 ]

1 голос
/ 13 октября 2011

Мой Matlab довольно ржавый, но если я вас правильно понял, вы, вероятно, захотите передать последовательность 1:Ne или 1:4000 для аргумента x в lowess, так как это должно быть x координаты точек, которые вы сглаживаете. Это предполагает, что вы предполагаете, что ваши очки действительно распределены равномерно.

Примерно так будет работать:

#Example matrix
M <- matrix(runif(1600),40,40)

#Smooth rows; transpose when smoothing over rows
M1 <- t(apply(M,1,FUN = function(x){lowess(1:length(x),x)$y}))

#Smooth columns; but don't transpose; fills by column already
M2 <- apply(M1,2,FUN = function(x){lowess(1:length(x),x)$y})

Я не включал разные промежутки, но вы можете добавить эти детали в себя.

Редактировать

Ах, но если вы ищете скорость, вам, вероятно, следует использовать loess.smooth напрямую. loess использует интерфейс формул, поэтому вам нужно вызвать loess.smooth напрямую. По умолчанию он отличается от lowess, поэтому будьте осторожны. Поменяв эту функцию, я сократил время работы почти на 1/4.

0 голосов
/ 13 октября 2011

Что ж, хотя ответ, предоставленный Джораном, не сработал так, как предполагалось, и Джоран потратил довольно много времени, пытаясь его отладить, что я ценю, в конечном итоге я не смог заставить это решение работать, и я недействительно знаю почему, но он не давал результатов, сравнимых с кодом Matlab.После возни я нашел это решение (которое может быть неоптимальным, но работает как положено)

loesscontrol=loess.control(surface="interpolate", statistics="approximate", trace.hat="exact", iterations=1)
spanNe=100/Ne
spanT=50/nsteps
spanNi=30/Ni

for (i in 1:nsteps){
        x<-1:Ne
        y<-spike_sum[1:Ne,i]
        smoothingNe<-loess.smooth(x, y, span=spanNe, degree=1, family="gaussian", evaluation=Ne)
        smooth_sumNe[1:Ne,i]<-smoothingNe$y
}

for (i in 1:Ne){
        x<-1:nsteps
        y<-smooth_sumNe[i,1:nsteps]
        smoothingT<-loess.smooth(x, y, span=spanT, degree=1, control=loesscontrol, family="gaussian", evaluation=nsteps)
        smoother_sumNe[i,1:nsteps]<-smoothingT$y        
}

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

...