配件定制分布时避免数值问题

指定适当的初始参数值

因为审查数据的值并不清楚,最大似然估计需要更多的信息。特别是,概率密度函数(pdf),累积分布函数(cdf),和足够的初始参数值是计算loglikelihood所必需的。您可以使用evpdf和evcdf函数来指定pdf和运作。

首先,生成一些未经审查的极端值数据。

rng (0,“旋风”);n = 50;μ= 5;σ= 2.5;x = evrnd(μ、σ,n, 1);

接下来,审查的任何值大于预定的截止代替这些值截断值。这种类型的审查被称为II型审查。

c = (x > 7);x (c) = 7;

检查的比例审查观察。

和(c) / (c)长度

ans = 0.1200

百分之十二的原始数据是right-censored截止7点。

画一个柱状图的数据,包括堆叠柱形代表审查观察。

[uncensCnts,边缘]= histcounts (x (~ c), 10);censCnts = histcounts (x (c),边);栏((1:end-1) + diff边缘(边缘)/ 2,[uncensCnts ' censCnts '],“堆叠”)传说(“完全观测数据”,审查数据的,“位置”,“西北”)

尽管观测数据包括审查,审查的部分观测相对较小。因此,时刻的方法可以为参数估计提供合理的起点。计算的初始参数值μ和σ对应于观测到未经审查的数据的平均值和标准偏差。

sigma0 = sqrt(6) *性病(x (~ c)) /π

sigma0 = 2.3495

mu0 =意味着(x ~ (c))ψ(1)* sigma0

mu0 = 3.5629

找到两个极值分布的极大似然估计参数,以及近似95%置信区间。指定审查向量,pdf, cdf实验组,初始参数值。因为σ(尺度参数)必须是积极的,您还需要指定参数范围较低。

[paramEsts, paramCIs] =大中型企业(x,“审查”c…“pdf”@evpdf,“提供”@evcdf,…“开始”,(mu0 sigma0),下界的(从0))

paramEsts =1×24.5530 - 3.0215

paramCIs =2×23.6455 2.2937 5.4605 3.7494

指定logpdf和logsf

安装一个定制的初始分配需要一个参数,并确定如何好或坏的起点是一个先天是很困难的。如果您指定一个起点,远离最大似然估计,一些观察可以位于远离的反面极端值分布对应的起点。在这种情况下,会发生这些情况之一:

提供价值也可能成为如此之小,它下溢,但这种情况并不构成问题。

条件可以导致问题大中型企业计算loglikelihood,因为每个导致loglikelihood值负的优化算法大中型企业不能处理。

检查发生了什么有不同的起点。

开始= [1];试一试[paramEsts, paramCIs] =大中型企业(x,“审查”c…“pdf”@evpdf,“提供”@evcdf,…“开始”开始,下界的(从0))抓我disp (ME.message)结束

定制的累积分布函数返回值大于或等于1。

在这种情况下,第二个问题状况发生。一些提供值的初始参数估计是1,所以loglikelihood是无限的。你可以尝试设置FunValCheck控制参数从通过使用选项名称-值参数。的从选择禁用检查非限定的似然值。然而,解决这个数值问题的最好方法是在其根。

极值cdf的形式

p = 1 - exp (exp ((xμ)。/σ))

loglikelihood审查观测的贡献是他们的生存函数的日志(SF)值,或日志(1-cdf)。极端值分布,科幻的日志exp ((xμ)。/σ)。如果你计算loglikelihood直接使用日志科幻小说,而不是计算日志(1 - (1-exp (logSF))),你就可以避免cdf实验组的舍入的问题。观察提供的值是不区分1双精度的日志科幻值,很容易能上演的非零值。例如,一个提供的价值(1-1e-20)轮,1在双精度,因为双精度每股收益是关于2 e-16。

SFval = 1 e-20;cdfval = 1 - SFval

cdfval = 1

该软件可以很容易地代表相应的日志科幻小说的价值。

日志(SFval)

ans = -46.0517

同样的情况也适用于日志pdf格式;未经审查的观测到的贡献loglikelihood pdf的记录值。您可以使用日志pdf直接,而不是计算日志(exp (logpdf))为了避免下溢问题,pdf不区分从零双精度。该软件可以很容易地将日志pdf表示为一个有限的负数。例如,一个pdf的价值1 e - 400在双精度下溢,因为双精度最小正浮点数是关于2 e - 308。

logpdfval = -921;pdfval = exp (logpdfval)

pdfval = 0

使用大中型企业功能,您可以指定一个自定义分布与日志pdf和日志科幻(而不是pdf和cdf)通过设置logpdf和logsf名称-值参数。与pdf和提供功能,日志pdf和日志科幻没有内置函数。因此,您需要创建匿名函数,计算这些值。

evlogpdf = @ (x,μ、σ)((xμ)。/σ- exp ((xμ)。/σ))——日志(σ);evlogsf = @ (x,μ、σ)exp ((xμ)。/σ);

使用相同的起点,备用日志pdf和日志科幻的极端值分布规范让问题解决。

开始= [1];[paramEsts, paramCIs] =大中型企业(x,“审查”c…“logpdf”evlogpdf,“logsf”evlogsf,…“开始”开始,下界的(从0))

paramEsts =1×24.5530 - 3.0215

paramCIs =2×23.6455 2.2937 5.4605 3.7494

这个过程并不总是解决问题的起点,所以建议仔细选择的起点。

供应梯度优化功能fmincon

默认情况下,大中型企业使用函数fminsearch找到参数值最大化的loglikelihood定制的分布。fminsearch导数免费使用一个优化算法,使它的一个不错选择这种类型的问题。然而,对于一些问题,选择一个优化算法,它使用loglikelihood函数的导数可以收敛于最大似然估计的区别,特别是在出发点是远离最终答案。提供衍生品也可以加快收敛。

您可以指定OptimFun名称-值参数在大中型企业作为fmincon使用fmincon函数(需要优化工具箱)。的fmincon功能包括优化算法,可以利用导数信息。利用算法fmincon,指定一个自定义使用loglikelihood分布函数,不仅返回loglikelihood写的,但是它的梯度。loglikelihood函数的梯度向量的偏导数的参数。

这个策略需要额外准备编写代码,同时计算loglikelihood及其梯度。定义一个函数的名字helper_evnegloglike在一个单独的文件。

函数[附近,ngrad] = helper_evnegloglike(参数x,岑、频率)% HELPER_EVNEGLOGLIKE负对数似极值%分布。%这个函数只支持时避免数值问题的金宝app例子%安装自定义分布(customdist2demo.mlx)和可能的变化%将来的版本。如果元素个数(params) ~ = 2错误消息(“统计数据:probdists: WrongParameterLength”2));结束μ= params (1);σ= params (2);nunc =总和(1-cens);z = (x -μ)。/σ;expz = exp (z);nll =和(expz)和(z(~岑))+ nunc。*日志(σ);如果nargout > 1 ngrad =[总和(expz)。/σ+ nunc. /σ,…总和(z。* expz)。/σ+ (z(~岑))。/σ+ nunc. /σ);结束

这个函数helper_evnegloglike回报的负loglikelihood因为价值观和梯度值大中型企业最大限度地减少负面loglikelihood。

计算基于最大似然估计使用梯度优化算法,指定nloglf,OptimFun,选项名称-值参数。nloglf指定自定义函数,计算负loglikelihood,OptimFun指定fmincon优化功能,选项指定fmincon使用第二个输出自定义函数的梯度。

开始= [1];[paramEsts, paramCIs] =大中型企业(x,“审查”c“nloglf”@helper_evnegloglike,…“开始”开始,下界的(负无穷,0),…“OptimFun”,“fmincon”,“选项”statset (“GradObj”,“上”))

paramEsts =1×24.5530 - 3.0215

paramCIs =2×23.6455 2.2937 5.4605 3.7493