r - 如何对使用 R 包 CausalImpact 中实现的贝叶斯结构时间序列模型进行的预测施加限制？

Question

CausalImpact来自 Google 的 R 包为 Bayesian structural time series modeling 提供了一个用户友好的界面方法，它允许人们在没有真实控制组的情况下获得反事实预测并估计因果影响(例如，广告事件的效率)。我在使用该软件时面临的问题是，在某些情况下，当建模响应是计数变量时，预测区间和/或预测平均值可能会呈现负值，这对于此类变量显然是不可信的。

一种直接的解决方法是使用响应的对数转换值，然后将结果转换回原始比例(实际上是包作者的 this approach has already been mentioned)。然而，在解释原始尺度上的结果时，转换回由 CausalImpact 产生的任何汇总统计数据并不是很有帮助。这就是我的意思 - 考虑以下示例:

y = c(7, 18, 11, 3, 3, 2, 89, 94, 48, 74,
      21, 13, 5, 9, 10, 18, 12, 4, 8, 4, 12, 8, 6, 7, 6)
x = c(7, 22, 28, 13, 16, 6, 4, 2, 2, 24, 8,
      9, 5, 5, 8, 7, 5, 11, 3, 4, 5, 1, 4, 2, 6)

dat = as.ts(cbind(y, x))
pre.period <- c(1, 6)
post.period <- c(7, 25)

library(CausalImpact)
impact <- CausalImpact(dat, pre.period, post.period)
plot(impact)

summary(impact)
Posterior inference {CausalImpact}

                         Average        Cumulative     
Actual                   24             448            
Prediction (s.d.)        2.9 (2.1)      54.8 (39.1)    
95% CI                   [-0.98, 6.8]   [-18.69, 129.1]

Absolute effect (s.d.)   21 (2.1)       393 (39.1)     
95% CI                   [17, 25]       [319, 467]     

Relative effect (s.d.)   718% (71%)     718% (71%)     
95% CI                   [582%, 852%]   [582%, 852%]   

Posterior tail-area probability p:   0.00111
Posterior prob. of a causal effect:  99.88901%

For more details, type: summary(impact, "report")

从上图和汇总表可以看出，Actual 的 95% 可信区间包含负值。干预的绝对效果为 21，累积效果为 393，相对效果为 718%。

现在，让我们在对数转换的响应上重新运行相同的模型:

ylog = log(y)
dat2 = as.ts(cbind(ylog, x))
impactLog <- CausalImpact(dat2, pre.period, post.period)
plot(impactLog)

summary(impactLog)
Posterior inference {CausalImpact}

                         Average       Cumulative   
Actual                   2.6           49.2         
Prediction (s.d.)        1.1 (0.28)    20.1 (5.29)  
95% CI                   [0.54, 1.6]   [10.17, 30.3]

Absolute effect (s.d.)   1.5 (0.28)    29.2 (5.29)  
95% CI                   [1, 2.1]      [19, 39.1]   

Relative effect (s.d.)   145% (26%)    145% (26%)   
95% CI                   [94%, 194%]   [94%, 194%]  

Posterior tail-area probability p:   0.00111
Posterior prob. of a causal effect:  99.88901%

For more details, type: summary(impact, "report")

当然，试图解释上表中反对数转换的影响相关汇总统计数据在原始加法量表上是没有意义的。例如，绝对效应变为 exp(1.5) = 4.48，累积效应变为 exp(49.2) = 2.33e+21，等等。反对数变换必须对用于构建上述汇总统计的原始值进行，但是我似乎不知道这些原始值在哪里，以及如何在原始尺度上计算有意义的汇总统计数据。任何帮助将不胜感激。

Answer 1

对于效果，您可以进行反对数变换 response和 point.pred然后减去以获得效果:

您的原始非日志结果:

(平均，累计)

> impact$summary$AbsEffect
[1]  20.60043 391.40809

要使用对数转换数据获得平均效果:

> mean(exp(impactLog$series$response[7:25,]) - exp(impactLog$series$point.pred[7:25,]))
[1] 20.44092

...以及累积效应:

> sum(exp(impactLog$series$response[7:25,]) - exp(impactLog$series$point.pred[7:25,]))
[1] 388.3774

(两个结果之间的微小差异是由于 MCMC 过程中的随机性。)

如果你想要逐点效果，它和上面一样，但没有 mean或 sum .

对于绝对平均效应的可信区间:

> # lower
> mean(exp(impactLog$series$response[7:25,]) - exp(impactLog$series$point.pred.upper[7:25,]))
[1] 15.22143

> # upper
> mean(exp(impactLog$series$response[7:25,]) - exp(impactLog$series$point.pred.lower[7:25,]))
[1] 22.40529

对于累积可信区间，替换 mean与 sum .