预测平均值和标准差答案

【问题标题】：Forecast Mean and Standard Deviation预测平均值和标准差
【发布时间】：2014-03-21 18:07:10
【问题描述】：

如果这是一个简单的问题，我深表歉意，但在过去的一周里我无法找到任何答案，这让我发疯了。

背景信息：我有一个数据集，用于跟踪 5 个人在 5 年内的体重。每年，我都有一个群体中个体的体重分布，我从中计算平均值和标准差。数据如下：

 Year = [2002,2003,2004,2005,2006]
 Weights_2002 = [12, 14, 16, 18, 20]
 Weights_2003 = [14, 16, 18, 20,20]
 Weights_2004 = [16, 18, 20, 22, 18]
 Weights_2005 = [18, 21, 22, 22, 20]
 Weights_2006 = [2, 21, 19, 20, 20]

问题：我如何预测该组未来 10 年的年度体重分布？理想情况下，我希望随着时间的推移，均值的不确定性会增加。同样，我希望标准偏差的不确定性也增加。换一种说法，我想预测未来的权重分布，兼顾两者：

数据中的自然方差
不确定性增加。

任何帮助将不胜感激。如果有人可以建议如何在 R 中执行此操作，那就更好了。

谢谢大家！

【问题讨论】：

这两个资源statmethods.net/advstats/timeseries.html 和robjhyndman.com/software/forecast 提供了有关使用 R 进行预测的出色教程。
@hrbrmstr 我之前都看过，但似乎无法找到可以解释数据自然变化和增加不确定性的预测。也许我遗漏了一些东西......我还使用 R 中的 Holt-Winters 包进行了投影，但我相信只能输入数据的平均值。它似乎没有考虑自然差异。但感谢您的快速响应和建议
您想要基于整个组的预测还是针对每个人的预测？
我想要整个团队的预测。我想预测该组的权重分布将如何随时间变化。就像我有 2006 年该组的体重分布一样，我也想预测 2007 年该组的体重分布。希望这很清楚
这似乎是分层分析的应用程序，您可以在其中对组趋势和方差进行建模。但我不确定您如何将未来的“不确定性”建模为时间的函数。我想这些信息不包含在现有数据中，因此必须从外部提供。

标签： r statistics distribution regression forecasting

【解决方案1】：

没有关于如何使用 R 中的预测工具的具体建议，即您的问题的 cmets，这是使用蒙特卡洛模拟的替代方法。

首先，一些内务处理：Weights_2006 中的值 2 要么是拼写错误，要么是异常值。由于我不知道是哪一个，我会假设它是一个异常值并将其排除在分析之外。

其次，你说你想基于increasing uncertainty 来预测分布。但您的数据不支持这一点。

Year <- c(2002,2003,2004,2005,2006)
W2   <- c(12, 14, 16, 18, 20)
W3   <- c(14, 16, 18, 20,20)
W4   <- c(16, 18, 20, 22, 18)
W5   <- c(18, 21, 22, 22, 20)
W6   <- c(NA, 21, 19, 20, 20)
df   <- rbind(W2,W3,W4,W5,W6)
df   <- data.frame(Year,df)

library(reshape2)   # for melt(...)
library(ggplot2)
data <- melt(df,id="Year", variable.name="Individual",value.name="Weight")
ggplot(data)+
  geom_histogram(aes(x=Weight),binwidth=1,fill="lightgreen",colour="grey50")+
  facet_grid(Year~.)

平均权重随着时间的推移而增加，但方差减少。查看各个时间序列即可了解原因。

ggplot(data, aes(x=Year, y=Weight, color=Individual))+geom_line()

一般来说，一个人的体重会随着时间线性增加（大约每年 2 个单位），直到达到 20 时，它会停止增加但会波动。由于您的初始分布是均匀的，因此体重较轻的人会随着时间的推移而增加，从而推动平均值上升。但是较重的人的体重停止增长。所以分布在 20 左右“聚集”起来，导致方差减小。我们可以从数字中看出这一点：均值增加，标准差减小。

smry <- function(x)c(mean=mean(x),sd=sd(x))
aggregate(Weight~Year,data,smry)
#   Year Weight.mean  Weight.sd
# 1 2002  16.0000000  3.1622777
# 2 2003  17.6000000  2.6076810
# 3 2004  18.8000000  2.2803509
# 4 2005  20.6000000  1.6733201
# 5 2006  20.0000000  0.8164966

我们可以使用蒙特卡罗模拟来模拟这种行为。

set.seed(1)
start <- runif(1000,12,20)
X <- start
result <- X
for (i in 2003:2008){
  X <- X + 2
  X <- ifelse(X<20,X,20) +rnorm(length(X))
  result <- rbind(result,X)
}
result <- data.frame(Year=2002:2008,result)

在这个模型中，我们从 1000 个人开始，他们的体重在 12 到 20 之间形成均匀分布，就像您的数据一样。在每个时间步，我们将权重增加 2 个单位。如果结果大于 20，我们将其裁剪为 20。然后我们添加分布为 N[0,1] 的随机噪声。现在我们可以绘制分布图了。

model <- melt(result,id="Year",variable.name="Individual",value.name="Weight")
ggplot(model,aes(x=Weight))+
  geom_histogram(aes(y=..density..),fill="lightgreen",colour="grey50",bins=20)+
  stat_density(geom="line",colour="blue")+
  geom_vline(data=aggregate(Weight~Year,model,mean), aes(xintercept=Weight), colour="red", size=2, linetype=2)+
  facet_grid(Year~.,scales="free")

红色条显示每年的平均体重。

如果您认为个体体重的自然变化会随着时间的推移而增加，则使用N[0,sigma] 作为模型中的误差项，sigma 随着Year 的增加而增加。问题是您的数据中没有任何内容可以支持这一点。

【讨论】：