least absolute deviation matlab

dav
2013 7 2
1 个回答
回答已采纳
9 次查看（30 天）

0 个投票

Hi,
Is the a way to do parameter estimation using Least Absolute deviation with constraints in matlab.
If so, please let me know of a place to read about it.
thanks.

请先登录，再回答此问题。

 采纳的回答

Matt J
Matt J 2013-7-2

1 个投票

If you are talking about the problem
min_x sum( abs(r(x)))
it is equivalent to
min_{x,d} sum(d)
s.t.
r(x) <=d
-r(x)<=d
The above formulation is differentiable if r(x) is differentiable, so FMINCON can handle it. You can obviously add any additional constraints you wish, so long as they too are differentiable.

14 个评论
显示 12更早的评论 隐藏 12更早的评论

dav
dav 2013-7-2
编辑：dav 2013-7-2
Hi,
Thanks for your answer.
I have a dataset yt (all positive - size (100,1)). I need to fit a straight line to this data with the constraints that intercept parameter is greater than zero and slope parameter is between 0 1nd 1.
Could you please help me with this code? It is a little difficult to figure it out using your answer.
Thanks much!
Matt J
Matt J 2013-7-2
编辑：Matt J 2013-7-2
So, the problem is as follows?
min.
f(m,b) = sum( abs(yi-m*xi-b))
s.t.
b>=0
0<=m<=1
dav
dav 2013-7-2
编辑：dav 2013-7-2
yes, but
b > 0
0<= m < 1
Thanks.
Matt J
Matt J 2013-7-2
编辑：Matt J 2013-7-2
You cannot impose strict inequality constraints. The minimum can be ill-defined, for example in the minimization of f(a)=a^2 over a>0. If a=0 is not allowed, then what do you consider the minimizing point?
dav
dav 2013-7-2
编辑：dav 2013-7-2
Thanks. can I have something like a >0.00000001. Would you please help me with the code ?
dav
dav 2013-7-2
I can provide the code which generates the data. but I dont know how to estimate it?
dav
dav 2013-7-2
The data set is yt.
clc;
clear;
T = 300;
a0 = 0.1; a1 = 0.4;
ra = zeros(T+2000,1);
seed=123;
rng(seed);
ra = trnd(5,T+2000,1);
ytn = [];
epsi=zeros(T+2000,1);
simsig=zeros(T+2000,1);
unvar = a0/(1-a1);
for i = 1:T+2000
if (i==1)
simsig(i) = unvar;
s=(simsig(i))^0.5;
epsi(i) = ra(i) * s;
else
simsig(i) = a0+ a1*(epsi(i-1))^2;
s=(simsig(i))^0.5;
epsi(i) = ra(i)* s;
end
end
epsi2 = epsi.^2;
yt = epsi2(2001:T+2000);
Matt J
Matt J 2013-7-3
编辑：Matt J 2013-7-15
I don't understand your code, but here's an example for you to study
%%Fake data
xi=0:10;
yi=2*xi+bt+randn(size(xi))/2;
%%Linprog parameters
N=length(xi);
e=ones(N,1);
f=[0,0,e.'];
A=[xi(:),e,-speye(N);-xi(:), -e, -speye(N)];
b=[yi(:);-yi(:)];
lb=zeros(N+2,1);
ub=inf(N+2,1); ub(1)=1;
%%Perform fit and test
p=linprog(f,A,b,[],[],lb,ub);
slope=p(1),
intercept=p(2),
yf=slope*xi+intercept;
plot(xi,yi,xi,yf)
dav
dav 2013-7-15
编辑：dav 2013-7-15
Thank you very much
Could you please explain the following codes
f=[0,0,e.'];
A=[xi(:),e,-speye(N);-xi(:), -e, -speye(N)];
b=[yi(:);-yi(:)];
and also bt in yi=2*xi+bt+randn(size(xi))/2;
thanks
Matt J
Matt J 2013-7-15
编辑：Matt J 2013-7-15
Recall that the reformulation of the problem requires additional constraints
r(x) <=d
-r(x)<=d
where r(x) is your residual. The A,b pair that you've cited is the implementation of these additional constraints. bt is the true intercept of the simulated line data, yi. You can choose any value for it that you want, for simulation purposes, and see if the p(2) returned by linprog ends up being a good estimate of it.
dav
dav 2013-7-18
is it possible to add the constraint that the parameter estimates should be positive here?
Matt J
Matt J 2013-7-18
Not sure what you mean by "add". I already put positivity constraints in my example using the lb input argument to linprog. Similarly, I used ub to impose the upper bounds you originally mentioned.
dav
dav 2013-7-18
编辑：dav 2013-7-18
thank you.
I used your code to test the data set I have. BOTH parameter estimates should be 0.1. However, the estimates I get using your code are very different. is there a way to fix it. From a standard theory I know that when you regress the y vector I have mentioned on the x vector I should get parameters, both equal to 0.1
clc;
clear;
p=1;
T = 300;
a0 = 0.1; a1 = 0.1;
seed=123;
ra = randn(T+2000,1);
epsi=zeros(T+2000,1);
simsig=zeros(T+2000,1);
unvar = a0/(1-a1);
for i = 1:T+2000
if (i==1)
simsig(i) = unvar;
s=(simsig(i))^0.5;
epsi(i) = ra(i) * s;
else
simsig(i) = a0+ a1*(epsi(i-1))^2;
s=(simsig(i))^0.5;
epsi(i) = ra(i)* s;
end
end
epsi2 = epsi.^2;
y = epsi2(2001:T+2000); % THIS IS THE DATA SET I WANT TO TEST.
len = length(y);
x = zeros(len,p);
for i = 1:p
x(1+i:len,i) = y(1:len-i,1); % THIS IS THE x VECTOR
end
N=length(x);
e=ones(N,1);
f=[0,0,e.'];
A=[x(:),e,-speye(N);-x(:), -e, -speye(N)];
b=[y(:);-y(:)];
lb=zeros(N+2,1);
ub=inf(N+2,1); ub(1)=1;
%%Perform fit and test
p=linprog(f,A,b,[],[],lb,ub);
slope=p(1),
intercept=p(2),
Matt J
Matt J 2013-7-18
编辑：Matt J 2013-7-18
All I can say is that I tested the example code I gave you. It worked fine for me and produced accurate estimates.

请先登录，再进行评论。

更多回答（0 个）

类别

在 帮助中心 和 File Exchange 中查找有关 Linear Programming and Mixed-Integer Linear Programming 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by