为什么要使用极值理论?

机器算法验证 分位数 极值
2022-01-26 01:15:06

我来自土木工程,我们使用极值理论,例如GEV分布来预测某些事件的值,例如最大风速,即98.5%的风速会低于该值。

我的问题是为什么要使用这样的极值分布如果我们只使用整体分布并获得98.5% 概率的值,会不会更容易

4个回答

免责声明:在以下几点,这大致假定您的数据是正态分布的。如果您实际上是在设计任何东西,那么请与强大的统计专家交谈,并让该人在线上签名,说明级别将是多少。与其中 5 人或 25 人交谈。此答案适用于询问“为什么”的土木工程专业学生,而不适用于询问“如何”的工程专业人士。

我认为问题背后的问题是“什么是极值分布?”。是的,它是一些代数 - 符号。所以呢?对?

让我们想想 1000 年的洪水。他们很大。

当它们发生时,它们会杀死很多人。许多桥梁正在倒塌。
你知道哪座桥不会倒塌吗?我愿意。你还没有……还没有。

问题:哪座桥不会在 1000 年的洪水中倒塌?
答:设计的桥梁可以承受它。

您需要按照自己的方式进行操作的数据:
因此,假设您拥有 200 年的每日水数据。那里有1000年的洪水吗?不是远程。您有一个分布尾部的样本。你没有人口。如果您了解所有洪水的历史,那么您将拥有全部数据。让我们考虑一下。您需要多少年的数据,多少个样本才能获得至少一个可能性为千分之一的值?在一个完美的世界中,您至少需要 1000 个样本。现实世界是混乱的,所以你需要更多。您开始在大约 4000 个样本处获得 50/50 的赔率。您开始得到保证在大约 20,000 个样本中拥有超过 1 个样本。样本并不意味着“一秒钟与下一秒钟的水”,而是对每个独特变异来源的衡量——比如逐年变化。一项措施超过一年,与另一年的另一项措施一起构成两个样本。如果您没有 4,000 年的良好数据,那么您可能在数据中没有 1000 年洪水的示例。好消息是 - 你不需要那么多数据来获得好的结果。

以下是如何用更少的数据获得更好的结果:
如果您查看年度最大值,您可以将“极值分布”拟合到 year-max-levels 的 200 个值,您将得到包含 1000 年洪水的分布-等级。它将是代数,而不是实际的“它有多大”。您可以使用该等式来确定 1000 年洪水的规模。然后,鉴于水的体积 - 你可以建造你的桥梁来抵抗它。不要追求精确的值,而要追求更大的值,否则您将其设计为在 1000 年的洪水中失败。如果您是大胆的,那么您可以使用重新采样来确定您需要将其构建到的确切 1000 年值超出多少才能使其抵抗。

这就是为什么 EV/GEV 是相关分析形式的原因:
广义极值分布是关于最大值变化的程度。最大值的变化与平均值的变化完全不同。正态分布通过中心极限定理描述了许多“中心趋势”。

程序:

  1. 做以下 1000 次:
    i. 从标准正态分布中挑选 1000 个数字
    ii. 计算该组样本的最大值并将其存储
  2. 现在绘制结果的分布

    #libraries
    library(ggplot2)
    
    #parameters and pre-declarations
    nrolls <- 1000
    ntimes <- 10000
    store <- vector(length=ntimes)
    
    #main loop
    for (i in 1:ntimes){
    
         #get samples
         y <- rnorm(nrolls,mean=0,sd=1)
    
         #store max
         store[i] <- max(y)
    }
    
    #plot
    ggplot(data=data.frame(store), aes(store)) + 
         geom_histogram(aes(y = ..density..),
                        col="red", 
                        fill="green", 
                        alpha = .2) + 
         geom_density(col=2) + 
         labs(title="Histogram for Max") +
         labs(x="Max", y="Count")
    

这不是“标准正态分布”: 在此处输入图像描述

峰值为 3.2,但最大值上升至 5.0。它有偏斜。它不会低于大约 2.5。如果您有实际数据(标准法线)并且您只选择尾部,那么您将沿着这条曲线均匀随机地选择一些东西。如果你幸运的话,你会朝着中心而不是下尾巴。工程与运气相反——它是关于每次都能始终如一地实现预期的结果。随机数太重要了,不能任凭运气”(见脚注),尤其是对于工程师而言。最适合该数据的分析函数族 - 分布的极值族。

样本拟合:
假设我们有 200 个来自标准正态分布的年度最大值的随机值,我们将假装它们是我们 200 年的最高水位历史(无论这意味着什么)。要获得分布,我们将执行以下操作:

  1. 对“存储”变量进行采样(以制作简短/简单的代码)
  2. 拟合广义极值分布
  3. 找到分布的平均值
  4. 使用 bootstrapping 来找到均值变化的 95% CI 上限,因此我们可以针对我们的工程进行定位。

(代码假定上面已经先运行)

library(SpatialExtremes) #if it isn't here install it, it is the ev library
y2 <- sample(store,size=200,replace=FALSE)  #this is our data

myfit <- gevmle(y2)

这给出了结果:

> gevmle(y2)    
       loc      scale      shape     
 3.0965530  0.2957722 -0.1139021     

这些可以插入到生成函数中以创建 20,000 个样本

y3 <- rgev(20000,loc=myfit[1],scale=myfit[2],shape=myfit[3])

建立以下条件将使任何一年的失败几率为 50/50:

平均值(y3)
3.23681

这是确定 1000 年“洪水”水平的代码:

p1000 <- qgev(1-(1/1000),loc=myfit[1],scale=myfit[2],shape=myfit[3])
p1000

建立在以下基础上应该会给你 50/50 的失败几率在 1000 年的洪水中失败。

p1000
4.510931

为了确定 95% 的 CI 上限,我使用了以下代码:

myloc <- 3.0965530
myscale <- 0.2957722
myshape <- -0.1139021

N <- 1000
m <- 200
p_1000 <- vector(length=N)
yd <- vector(length=m)

for (i in 1:N){

      #generate samples
    yd <- rgev(m,loc=myloc,scale=myscale,shape=myshape)

    #compute fit
    fit_d <- gevmle(yd)

    #compute quantile
    p_1000[i] <- qgev(1-(1/1000),loc=fit_d[1],scale=fit_d[2],shape=fit_d[3])

}

mytarget <- quantile(p_1000,probs=0.95)

结果是:

> mytarget
     95% 
4.812148

这意味着,为了抵御 1000 年的大部分洪水,鉴于您的数据非常正常(不太可能),您必须为...

> out <- pgev(4.812148,loc=fit_d[1],scale=fit_d[2],shape=fit_d[3])
> 1/(1-out)

或者

> 1/(1-out)
   shape 
1077.829 

... 1078 年洪水。

底线:

  • 你有一个数据样本,而不是实际的总人口。这意味着您的分位数是估计值,并且可能会关闭。
  • 像广义极值分布这样的分布被构建为使用样本来确定实际的尾部。与使用样本值相比,它们在估计方面的情况要好得多,即使您没有足够的样本用于经典方法。
  • 如果你很健壮,上限就很高,但结果是——你不会失败。

祝你好运

PS:

  • 我听说一些土木工程设计的目标是 98.5%。如果我们计算的是第 98.5 个百分位数而不是最大值,那么我们会发现一条具有不同参数的不同曲线。我认为这是为了建立一个 67 年的风暴。
    1/(10.985)67
    imo 那里的方法是找到 67 年风暴的分布,然后确定平均值附近的变化,并获得填充,以便将其设计为在第 67 年风暴中成功而不是在其中失败。
  • 鉴于前一点,平民平均每 67 年就应该重建一次。因此,考虑到土木结构的使用寿命(我不知道那是什么),以每 67 年的全部工程和建设成本计算,在某个时候,在更长的风暴间期进行工程设计可能会更便宜。可持续的民用基础设施旨在至少持续一个人的寿命而不会失败,对吗?

PS:更有趣 - 一个 youtube 视频(不是我的)
https://www.youtube.com/watch?v=EACkiMRT0pc

脚注:Coveyou, Robert R. “随机数生成太重要了,不能靠运气。” 应用概率和蒙特卡罗方法以及动力学的现代方面。应用数学研究 3(1969):70-111。

如果您只对尾巴感兴趣,那么将数据收集和分析工作集中在尾巴上是有意义的。这样做应该更有效率。我强调了数据收集,因为在提出 EVT 分布的论点时,这方面经常被忽略。事实上,收集相关数据来估计某些领域的整体分布可能是不可行的。我将在下面更详细地解释。

如果您像@EngrStudent 的示例一样查看千分之一的洪水,那么要构建正态分布的主体,您需要大量数据来填充观察结果。您可能需要过去数百年来发生的每一次洪水。

现在停下来想想到底什么是洪水?当我的后院在一场大雨后被洪水淹没时,是洪水吗?可能不是,但是将洪水与非洪水事件划定的界限到底在哪里?这个简单的问题突出了数据收集的问题。您如何确保我们按照相同的标准收集身体上的所有数据几十年甚至几个世纪?要收集洪水分布体的数据几乎是不可能的。

因此,这不仅仅是分析效率的问题,而是数据收集的可行性问题:是对整个分布建模还是仅对尾部建模?

自然,有了尾巴,数据收集就容易多了。如果我们为大洪水定义足够高的阈值,那么我们就有更大的机会以某种方式记录所有或几乎所有此类事件。很难错过一场毁灭性的洪水,如果有任何文明存在,就会保存一些关于该事件的记忆。因此,鉴于数据收集在极端事件上比在可靠性研究等许多领域中的非极端事件上更加稳健,因此构建专门关注尾部的分析工具是有意义的。

您使用极值理论从观察到的数据中进行推断。通常,您拥有的数据不足以为您提供尾部概率的合理估计。以@EngrStudent 的千分之一年事件为例:这对应于找到分布的 99.9% 分位数。但是,如果您只有 200 年的数据,则最多只能计算 99.5% 的经验分位数估计值。

极值理论让您通过对尾部分布的形状做出各种假设来估计 99.9% 的分位数:它是平滑的,它以某种模式衰减,等等。

您可能会认为 99.5% 和 99.9% 之间的差异很小;毕竟只有0.4%。但这是概率上的差异,当你处于尾部时,它可以转化为分位数的巨大差异。这是伽马分布的示意图,随着这些事情的发展,它没有很长的尾巴。蓝线对应 99.5% 分位数,红线对应 99.9% 分位数。虽然在垂直轴上它们之间的差异很小,但在水平轴上的分离却很大。对于真正的长尾分布,这种分离只会变得更大;伽马实际上是一个相当无害的案例。

在此处输入图像描述

通常,基础数据(例如,高斯风速)的分布是针对单个样本点的。第 98 个百分位会告诉您,对于任何随机选择的点,该值有 2% 的机会大于第 98 个百分位。

我不是土木工程师,但我想你想知道的不是任何一天的风速超过某个数字的可能性,而是最大可能阵风的分布,比如说,一年的过程。在这种情况下,如果每日阵风最大值呈指数分布,那么您想要的是 365 天内最大阵风的分布……这就是极值分布要解决的问题。