matlab如何确定 ncomp,如何使用MATLAB进行分位数回归?

function b = rq_fnm(X, y, p)

% Construct the dual problem of quantile regression

% Solve it with lp_fnm

%

%

[m n] = size(X);

u = ones(m, 1);

a = (1 - p) .* u;

b = -lp_fnm(X', -y', X' * a, u, a)';

function y = lp_fnm(A, c, b, u, x)

% Solve a linear program by the interior point method:

% min(c * u), s.t. A * x = b and 0 < x < u

% An initial feasible solution has to be provided as x

%

% Function lp_fnm of Daniel Morillo & Roger Koenker

% Translated from Ox to Matlab by Paul Eilers 1999

% Modified by Roger Koenker 2000--

% More changes by Paul Eilers 2004

% Set some constants

beta = 0.9995;

small = 1e-5;

max_it = 50;

[m n] = size(A);

% Generate inital feasible point

s = u - x;

y = (A' \  c')';

r = c - y * A;

r = r + 0.001 * (r == 0);    % PE 2004

z = r .* (r > 0);

w = z - r;

gap = c * x - y * b + w * u;

% Start iterations

it = 0;

while (gap) > small & it < max_it

it = it + 1;

%   Compute affine step

q = 1 ./ (z' ./ x + w' ./ s);

r = z - w;

Q = spdiags(sqrt(q), 0, n, n);

AQ = A * Q;          % PE 2004

rhs = Q * r';        % "

dy = (AQ' \ rhs)';   % "

dx = q .* (dy * A - r)';

ds = -dx;

dz = -z .* (1 + dx ./ x)';

dw = -w .* (1 + ds ./ s)';

% Compute maximum allowable step lengths

fx = bound(x, dx);

fs = bound(s, ds);

fw = bound(w, dw);

fz = bound(z, dz);

fp = min(fx, fs);

fd = min(fw, fz);

fp = min(min(beta * fp), 1);

fd = min(min(beta * fd), 1);

% If full step is feasible, take it. Otherwise modify it

if min(fp, fd) < 1

% Update mu

mu = z * x + w * s;

g = (z + fd * dz) * (x + fp * dx) + (w + fd * dw) * (s + fp * ds);

mu = mu * (g / mu) ^3 / ( 2* n);

% Compute modified step

dxdz = dx .* dz';

dsdw = ds .* dw';

xinv = 1 ./ x;

sinv = 1 ./ s;

xi = mu * (xinv - sinv);

rhs = rhs + Q * (dxdz - dsdw - xi);

dy = (AQ' \ rhs)';

dx = q .* (A' * dy' + xi - r' -dxdz + dsdw);

ds = -dx;

dz = mu * xinv' - z - xinv' .* z .* dx' - dxdz';

dw = mu * sinv' - w - sinv' .* w .* ds' - dsdw';

% Compute maximum allowable step lengths

fx = bound(x, dx);

fs = bound(s, ds);

fw = bound(w, dw);

fz = bound(z, dz);

fp = min(fx, fs);

fd = min(fw, fz);

fp = min(min(beta * fp), 1);

fd = min(min(beta * fd), 1);

end

% Take the step

x = x + fp * dx;

s = s + fp * ds;

y = y + fd * dy;

w = w + fd * dw;

z = z + fd * dz;

gap = c * x - y * b + w * u;

%disp(gap);

end

function b = bound(x, dx)

% Fill vector with allowed step lengths

% Support function for lp_fnm

b = 1e20 + 0 * x;

f = find(dx < 0);

b(f) = -x(f) ./ dx(f);

这个是大师Roger Koenker写的分位数回归包---matlab程序 内点(Frisch-Newton)算法。

A basic version of the interior point (Frisch-Newton) algorithm for quantile regression developed for the R quantreg package is also available for matlab. A C++ translation of the algorithm is also available from Ron Gallant's libcpp library. This algorithm is described in Koenker and Portnoy (Statistical Science, 1997).

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值