文章插图
， 其中
文章插图
代表所有的观测数据 。
3）预测值估计：得到
文章插图
之后 ， 我们从该分布中对（α ， θ）进行采样 ， 再代入状态空间方程(1)中对y进行预测 ， 得到
文章插图
， 其中
文章插图
表示时间点n之后y的预测值 。
文章插图
图3-2：展示某结构化时序数据及其背后各个状态转换的过程 。状态α 包含 Local trend:
文章插图
（局部趋势); local level:
文章插图
(局部趋势的均值) 以及协变量x,
文章插图
表示所有的观测数据;
文章插图
表示根据状态模型得到的预测数据 。
文章插图
分别表示
文章插图
的标准差这些参数均通过MCMC的方式得到估计 。
四、模型应用与代码实现以上我们给出了BSTS模型及MCMC方法的简要理论推导及结果输出 ， 核心目的就是对观测值y做出预测 。接下来我们将介绍如何在因果推断场景中应用BSTS模型 。
在对政策的效果评估上 ， 我们核心想要的是观测对象“反事实值” ， 例如“如果没有这个广告投放 ， 用户的浏览情况会怎样？”相较于传统的PSM或SCM方法 ， BSTS胜在其能够对于时间序列数据进行效果评估；同时利用贝叶斯估计输出反事实值y的预测 ， 并给出预测值的置信区间 ， 能一定程度上降低反事实值预测的波动性 ， 提升效应评估的准确性与稳定性 。
在实践应用上 ， 可以通过谷歌开源的CausalImpact包来实现BSTS模型 ， 在Python/ target=_blank class=infotextkey>Python和R中均可调用 ， 具体代码实现详见参考文献[7][8] 。
文章插图
图4-1：展示执行代码时的三个步骤：训练BSTS模型；反事实值预测；计算因果效应值 ， 包括效应值的点估计及置信区间 。
4.1 代码实现
下面通过模拟数据展示代码的具体命令
import tensorflow as tf
import tensorflow_probability as tfpfrom causalimpact import CausalImpact# 模型初始化 - 自定义时间序列数据：def plot_time_series_components(ci):component_dists = tfp.sts.decompose_by_component(ci.model, ci.observed_time_series, ci.model_samples)num_components = len(component_dists)mu, sig = ci.mu_sig if ci.mu_sig is not None else 0.0, 1.0for i, (component, component_dist) in enumerate(component_dists.items()):component_mean = component_dist.mean().numpy()component_stddev = component_dist.stddev().numpy()# 自定义观测方程以及真实值y：def plot_forecast_components(ci):component_forecasts = tfp.sts.decompose_forecast_by_component(ci.model, ci.posterior_dist, ci.model_samples)num_components = len(component_forecasts)mu, sig = ci.mu_sig if ci.mu_sig is not None else 0.0, 1.0for i, (component, component_dist) in enumerate(component_forecasts.items()):component_mean = component_dist.mean().numpy()component_stddev = component_dist.stddev().numpy()# 生成模拟数据 ， 包括一个实验组数据（有干预）以及两条对照组数据（无干预）observed_stddev, observed_initial = (tf.convert_to_tensor(value=https://www.isolves.com/it/cxkf/jiagou/2023-09-12/1, dtype=tf.float32),tf.convert_to_tensor(value=0., dtype=tf.float32))level_scale_prior = tfd.LogNormal(loc=tf.math.log(0.05 * observed_stddev), scale=1, name='level_scale_prior')# 设置先验分布initial_state_prior = tfd.MultivariateNormalDiag(loc=observed_initial[..., tf.newaxis], scale_diag=(tf.abs(observed_initial) + observed_stddev)[..., tf.newaxis], name='initial_level_prior')# 设置先验分布 ll_ssm = tfp.sts.LocalLevelStateSpaceModel(100, initial_state_prior=initial_state_prior,level_scale=level_scale_prior.sample())#训练时序模型ll_ssm_sample = np.squeeze(ll_ssm.sample().numpy())# 整合数据x0 = 100 * np.random.rand(100)# 对照组1x1 = 90 * np.random.rand(100)# 对照组2y = 1.2 * x0 + 0.9 * x1 + ll_ssm_sample#生成真实值yy[70:] += 10#设置干预点data = pd.DataFrame({'x0': x0, 'x1': x1, 'y': y}, columns=['y', 'x0', 'x1'])
推荐阅读
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
- 小李子狂撩37岁俄罗斯超模，聚会到凌晨离开，女方还曾迷倒C罗侃爷！
- 克里斯·埃文斯倾情演绎，英雄的诞生与正义的战斗，引发观众深思
- 网红rose被吐槽，模仿菲利斯假借孝心赚流量，网友：找茬博眼球
- 提拉米苏和慕斯的区别,哪个好吃 提拉米苏和慕斯的区别
- 小托马斯：科比的敬业NBA第一，他决不会像詹姆斯杜兰特选择抱团
- 粉丝超千万的网红rose被吐槽，处处模仿菲利斯，假借孝心博流量
- 格瓦斯是什么 格瓦斯是什么意思
- 美媒将职业生涯场均23分以上的29人划分成6档，詹姆斯阿杜同一档
- 陈佩斯豪宅被曝光，4个孙子和他长得非常相似
- 五步让你掌握Python数据结构
