Деконволюция с помощью R (пакет decon и deamer)

В ответ на ваш последний комментарий: ошибка сдвигает наблюдаемые значения. Сигнал, который вы хотите развернуть, находится где-то между 0 и ~ 0,3, я думаю. Вот некоторый код, использующий deamer:

library(actuar) # for rinvweibull
library(deamer)
set.seed(123)
RT <- rinvweibull(30000, shape=5.53861156, scale=488)/1000
RT <- RT[RT<1.5]
noise <- 1/rnorm(30000, mean=0.0023853421, sd=0.0004784688)/1000
noise <- noise[noise<1.5]

ST <- deamerSE(RT, errors=noise, from=0, to=0.3)
plot(ST)

Это то, что вы получите, используя непараметрическую деконволюцию (независимо от реализации, пакетов и т. д.). К вашему сведению, ваше отношение сигнал/шум чрезвычайно низкое... то, что вы наблюдаете, на самом деле почти просто шум. Это сильно влияет на оценку интересующей вас плотности, особенно при использовании непараметрических методов, таких как попытка найти иголку в стоге сена. Вам следует пересмотреть оценку плотности и попытаться получить только несколько интересующих величин...

Удачи, Жюльен Штирнеманн

В вашем сообщении есть несколько проблем. Во-первых: в задачах непараметрической деконволюции вы обычно не «знаете» распределение «у». Скорее у вас есть образец «y», который, как вы предполагаете, наблюдается с аддитивным шумом, «x» не наблюдается. Не делается никаких предположений относительно «у» или «х», а только «шум». Ваша презентация, кажется, подразумевает, что вы рассматриваете параметрическую проблему (для которой ни deamer, ни decon не помогают). Во-вторых: будьте осторожны, вы рассматриваете нецентрированный шум... с которым может справиться деамер, но не декон. Вот пример кода:

library(decon) # for DeconPdf
library(deamer) # for deamerKE

set.seed(12345)
shape=5; scale=1; mu=0; sd=0.2

x <- rweibull(5000, shape=shape, scale=scale)
noise <- rnorm(5000, mean=mu, sd=sd)
y=x+noise
curve(dweibull(x,shape,scale),lwd=2, from = 0, to = 2)

est <- deamerKE(y, noise.type="Gaussian", mu=mu, sigma=sd, from=0, to=2)
lines(est)

estDecon <- DeconPdf(y, sd, error="normal", fft=TRUE)
lines(estDecon, lty=2)

legend('topright', lty=c(1,1,2), lwd=c(2,1,1), 
    legend=c("true", "deamerKE", "DeconPdf"))

Как видно из графика, даже при центрированном шуме (в моем примере mu=0) оценка лучше с деамером: это из-за адаптивной оценки. Вы, вероятно, могли бы получить аналогичные результаты с помощью decon, но вам пришлось бы настроить параметр пропускной способности, используя функции, предусмотренные для этого в пакете. Что касается параметров, которые вы указали, преобразования Фурье чрезвычайно «плоские». Это очень затрудняет для любой универсальной реализации выбор подходящего параметра полосы пропускания (либо адаптивно, как в deamer, либо с использованием оценки, как в decon). Игра с параметром пропускной способности в deconPdf тоже не помогает, вероятно, из-за числовых ограничений. Ваша проблема потребует тонкой настройки кода функций deamer, чтобы можно было исследовать более крупные наборы моделей. Это также значительно увеличило бы время оценки. Не могли бы вы лучше рассмотреть возможность какой-либо перепараметризации вашей проблемы?

Бест, Жюльен Штирнеманн

После вашего второго сообщения: я не уверен, что полностью понимаю вашу проблему. Однако, насколько я понимаю, 2 варианта:

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

2) Если РТ наблюдается в одном и.и.д. время выборки и выборки в другой независимой выборке, тогда да, единственный выход — рассмотреть оценку плотности деконволюции. Однако, несмотря на то, что вы сделали некоторые параметрические предположения, используя методы подгонки, вы действительно не знаете плотности RT или DT. Рассматривая DT как шум, ваша проблема такова: RT = ST + шум со вспомогательной выборкой i.i.d. шум, заданный вашим образцом времени дискриминации. Вы хотите оценить плотность ST, которая не наблюдается. Единственный пакет, который может выполнять деконволюцию в этой ситуации, это deamer (насколько мне известно) с функцией deamerSE. Если я правильно изложил вашу проблему, вы должны посмотреть примеры в руководстве. Я бы также рекомендовал использовать необработанные данные без преобразований (по крайней мере, при первом анализе). Пример:

deamerSE(RT, errors=DT)

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

Бест, Жюльен Штирнеманн