How to perform mathematical operation of arrays inside loop?

1 次查看(过去 30 天)
I know it is a very basic question, but I have not been able to solve it. I have an array X(28*1). Using values of X for the first time, I would calculate E=B-A*X where E=(352*1), A(352*28) and B=(352*1). Then if the sum of E is larger than 10^-40, The follwoing operations are done:
C=diag(E*E');
W=diag(C)^-1;
X=(A'*W*A)^-1*A'*W*B
Then the updated value of X will be used for E=B-A*X operation. I understand whole thing is required to be inside a loop. But I can not find the exact process to write the code. If anyone has any idea, It would be much appriciated. I have attached the X for 1st iteration, A, and B as .mat file.
  6 个评论
Bruno Luong
Bruno Luong 2023-7-19
If you want to do robust linear fit with outliers, I recommend you do L1 fitting
minimize norm(A*X - B, 1)
rather than weighted 2-norm (more robust than l2, but still not very robust).
This requires linear programing linprog of optimization toolbox.
Matt J
Matt J 2023-7-19
编辑:Matt J 2023-7-19
I recommend you do L1 fitting...This requires linear programing linprog of optimization toolbox.
But, I'd like to add that there would be no need to implement it from scratch. You could instead use minL1lin() from this FEX download,

请先登录,再进行评论。

回答(3 个)

Matt J
Matt J 2023-7-18
编辑:Matt J 2023-7-18
for i=1:N
E=B-A*X;
if norm(E,1)<1e-40, continue; end
C=E.^2;
X=lscov(A,B,1./C);
end
  3 个评论
Matt J
Matt J 2023-7-20
编辑:Matt J 2023-7-20
@ANANTA BIJOY BHADRA I notice that the second column of your A matrix is zero, which iof course makes A ill-conditioned.
A=load('A').A_primary_A;
[l,h]=bounds(A(:,2))
l = 0
h = 0
cond(A)
ans = 2.5955e+16
cond(A(:, [1,3:end])) %skip second column
ans = 1.4142
Bruno Luong
Bruno Luong 2023-7-20
编辑:Bruno Luong 2023-7-20
Good catch. But I beleive the solution of "\" or lscov, both based on QR decomposition would not be affected by 0 column and remain accurate, despite the warning. The solution corresponds to 0-column is simply nil.
A=rand(10,3);
b=rand(10,1);
As=A; As(:,end+1)=0;
w = 0.1+rand(size(b));
format long
A\b
ans = 3×1
0.305726431573425 0.655213115202102 -0.126807101744095
As\b
Warning: Rank deficient, rank = 3, tol = 4.365185e-15.
ans = 4×1
0.305726431573425 0.655213115202102 -0.126807101744095 0
lscov(A,b,w)
ans = 3×1
0.365558371139105 0.511882753808717 -0.133776410641447
lscov(As,b,w)
Warning: A is rank deficient to within machine precision.
ans = 4×1
0.365558371139105 0.511882753808717 -0.133776410641447 0

请先登录,再进行评论。


Mrutyunjaya Hiremath
E = B-A*X;
while(sum(E(:))> 10^-40)
C=diag(E*E');
W=diag(C)^-1;
X=(A'*W*A)^-1*A'*W*B;
end

Voss
Voss 2023-7-18
E = B-A*X;
while sum(E,'all') > 1e-40
C = diag(E*E');
W = diag(C)^-1;
X = (A'*W*A)^-1*A'*W*B;
E = B-A*X;
end

类别

Help CenterFile Exchange 中查找有关 Programming 的更多信息

产品


版本

R2023a

Community Treasure Hunt

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

Start Hunting!

Translated by