指定足够的初始参数值

因为截尾数据的值不确切，所以最大似然估计需要更多的信息。特别是，概率密度函数（pdf）、累积分布函数（cdf）和足够的初始参数值是计算对数似然的必要条件。你可以使用evpdf和evcdf函数指定PDF和cdf。

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

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

接下来，通过用截止值替换任何大于预定截止值的值来审查这些值。这种类型的审查称为第二类审查。

c=（x>7）；x（c）=7；

检查审查意见的百分比。

和(c) / (c)长度

ans=0.1200

百分之十二的原始数据被右删，截取率为百分之七。

绘制数据的直方图，包括一个堆叠的条形图来代表被删的观测结果。

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

虽然数据包含截尾观测，但截尾观测的分数相对较小。因此，矩量法可以为参数估计提供合理的起点。计算μ和σ与未审查数据的观测平均值和标准偏差相对应。

sigma0=sqrt（6）*标准（x（~c））/pi

sigma0 = 2.3495

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

mu0 = 3.5629

找出两个极值分布参数的最大似然估计，以及近似的95%置信区间。指定审查矢量、pdf、cdf和初始参数值。因为σ(缩放参数)必须为正，还需要指定参数下限。

[paramEsts, paramCIs] =大中型企业(x,“审查”C．．.“pdf”，@evpdf，“cdf”@evcdf,．．.“开始”, (mu0 sigma0),下界的(从0))

paramEsts =1×24.5530 3.0215

paramCIs =2×23.6455 2.2937 5.4605 3.7494

指定logpdf和logsf

拟合自定义分布需要对参数进行初始猜测，而先验地确定起点的好坏可能很困难。如果指定的起点距离最大似然估计值较远，则某些观测值可能位于与起点相对应的极值分布尾部较远的位置。在这种情况下，可能出现以下情况之一：

cdf值也可能变得如此之小，以致于溢出，但这种情况并不构成问题。

这两种情况都可能导致问题大中型企业计算对数似然，因为每一个都导致的对数似然值负，其中的优化算法大中型企业不能处理。

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

Start = [1 1];试一试[paramEsts, paramCIs] =大中型企业(x,“审查”C．．.“pdf”，@evpdf，“cdf”@evcdf,．．.“开始”开始下界的(从0))接住我disp (ME.message)结束

自定义累积分布函数返回值大于等于1。

在这种情况下，会出现第二个问题条件。初始参数猜测处的一些cdf值恰好为1，因此对数似然为无穷大。你可以试着设置FunValCheck控制参数关通过使用选项名称-值参数。的关选项禁用对非有限似然值的检查。然而，解决这个数值问题的最佳方法是从根本上解决。

极值cdf有如下形式

P = 1 -exp(-exp((x-mu)./sigma))

截尾观测值对对数似然的贡献是其生存函数(SF)值的对数，或日志(1-cdf). 对于极值分布，SF的对数为exp ((xμ)。/σ)．如果你直接使用log SF来计算对数可能性，而不是计算日志(1 - (1-exp (logSF)))，您可以避免使用cdf的舍入问题。cdf值不能区分的观测值1在双精度中，具有易于表示为非零值的对数SF值。例如，cdf值为（1-1e-20）轮,1在双精度中，因为双精度EPS是关于2e-16．

SFval=1e-20；cdfval=1-SFval

cdfval=1

软件可以很容易地表示相应SF值的日志。

日志(SFval)

ans = -46.0517

同样的情况也适用于日志pdf;未经审查的观测结果对对数似然的贡献是其PDF值的对数。您可以直接使用日志pdf，而不是计算日志(exp (logpdf))，以避免pdf在双精度下无法与零区分的底流问题。该软件可以很容易地将日志pdf表示为有限的负数。例如，pdf值为1 e - 400下流在双精度，因为双精度最小正浮点数是关于2 e - 308．

logpdfval = -921;pdfval = exp (logpdfval)

pdfval = 0

使用大中型企业函数，您可以通过设置logpdf和logsf名称值参数。与pdf和cdf函数不同，log pdf和log SF没有内置函数。因此，您需要创建计算这些值的匿名函数。

Evlogpdf = @(x,mu,sigma) (x-mu)/∑- exp((x)./∑)- log(∑);Evlogsf = @(x,mu,sigma) -exp((x-mu)./sigma);

使用相同的起点，极值分布的备用日志pdf和日志SF规范使得问题可以解决。

开始=[1]；[参数，参数]=mle（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查找使自定义分布的对数似然性最大化的参数值。fminsearch使用无导数的优化算法，使其成为这类问题的一个很好的选择。然而，对于一些问题，选择一个使用对数似然函数导数的优化算法可能会决定是否收敛于最大似然估计，特别是当起点离最终答案很远的时候。提供衍生品也可以加快收敛速度。

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

这种策略需要额外的准备工作，以便编写计算对数似然及其梯度的代码。定义一个名为helper_evnegloglike在一个单独的文件中。

函数[nll，ngrad]=helper_evnegloglike（参数，x，cens，freq）%极值的HELPER_EVNEGLOGLIKE负对数似然%分布。%此函数仅支持示例，以避免出现数金宝app字问题%拟合自定义发行版(customdist2demo.mlx)并可能在%未来的版本。如果元素个数(params) ~ = 2错误消息('stats:probdists:ErrorParameterLength',2));结束μ= params (1);σ= params (2);nunc =总和(1-cens);Z = (x -) /;expz = exp (z);NLL = sum(expz) - sum(z(~cens)) + NLL .*log(sigma);如果ngrad = [-sum(expz). ngrad = [expz] . ngrad = [expz] . ngrad = [expz] . ngrad = [expz]。/σ+ nunc. /σ,．．.-总和（z.*expz）。/sigma+sum（z（~cens））/sigma+nunc./sigma]；结束

这个函数helper_evnegloglike返回对数似然值和梯度值的负值，因为大中型企业最小化负对数似然。

要使用基于梯度的优化算法计算最大似然估计，请指定nloglf，OptimFun,选项名称-值参数。nloglf指定计算负对数似然的自定义函数，OptimFun指定fmincon为优化函数，和选项指定fmincon使用自定义函数的第二个输出作为渐变。

开始=[1]；[参数，参数]=mle（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