r - 使用 R 进行反卷积(decon 和 deamer 包)

Question

我有一个以下形式的模型:y = x + 噪声。我知道“y”和噪声的分布，并希望得到“x”的分布。所以我尝试用 R 对分布进行反卷积。我发现了 2 个包(decon 和 deamer)，我认为这两种方法应该或多或少是相同的，但我不明白为什么用 DeconPdf 进行反卷积会给我一个类似于正态分布的东西，并且用 deamerKE 去卷积给我一个均匀的分布。这是一个示例代码:

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

set.seed(12345)
y <- rweibull(10000, shape=5.780094, scale=0.00204918)
noise <- rnorm(10000, mean=0.002385342, sd=0.0004784688)
sdnoise <- sd(noise)

est <- deamerKE(y, noise.type="Gaussian", 
                mean(noise), sigma=sdnoise)
plot(est)

estDecon <- DeconPdf(y, sdnoise, error="normal", fft=TRUE)
plot(estDecon)

---

编辑(回应Julien Stirnemann):

我不确定重新参数化。我的实际问题是:我有反应时间 (RT)，理论上可以描述为 f(RT) = g(辨别时间) + h(选择时间)，其中 f、g 和 h 可以是这些时间值的转换。我的数据集中有“RT”和“判别时间”值。我对选择时间或 h(选择时间)感兴趣。通过核密度估计，我发现威 bool 分布最适合 1/RT 值，而正态分布最适合 1/(判别时间)。这就是为什么我可以将我的问题写成 1/RT = 1/(判别时间) + h(选择时间) 或 y = x + 噪声(我认为噪声为 1/(判别时间))。通过模拟这些 react 时间，我得到了具有以下参数的以下分布:

y <- rweibull(10000, shape=5.780094, scale=0.00204918)
noise <- rnorm(10000, mean=0.002385342, sd=0.0004784688)

重新参数化是什么意思？使用不同的值，例如对于比例参数？

Answer 1

您的帖子中存在几个问题。第一:在非参数反卷积问题中，您通常“不知道”“y”的分布。相反，你有一个你假设观察到的带有加性噪声的“y”样本，“x”是未观察到的。没有对“y”或“x”做出任何假设，而只对“噪音”做出假设。您的演示文稿似乎暗示您正在考虑一个参数问题(对此 deamer 或 decon 都没有任何帮助)。第二:要小心，您正在考虑非中心噪声……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)，deamer 的估计值也更好:这是因为自适应估计。虽然您可能会使用 decon 获得类似的结果，但是您必须使用包中为此提供的函数来调整带宽参数。关于您提供的参数，傅立叶变换非常“平坦”。这使得任何通用的实现都很难选择合适的带宽参数(无论是在 deamer 中自适应还是在 decon 中使用估计)。在 deconPdf 中使用带宽参数也无济于事，可能是因为数值限制。您的问题需要对 deamer 函数的代码进行一些微调，以允许探索更大的模型集合。这也会显着增加估计时间。您是否应该考虑以某种方式重新参数化您的问题？

最好的，朱利安施蒂内曼