Solving a non linear equation using the Least squares method
2 views (last 30 days)
Show older comments
This concern MRI data and diffusion tensor imaging
I am trying to estimate D within each voxel using the data in data1.mat and i am not sure if i have done it correctly
Any help would appreciated
syms d_xx d_yy d_zz d_xy d_xz d_yz real
D = [d_xx, d_xy, d_xz; d_xy, d_yy, d_yz; d_xz, d_yz, d_zz]
syms x y z real
g = [x; y; z]
sympref('FloatingPointOutput',false);
syms b S_0
S = S_0*exp(-b*g'*D*g)
log(S)
logS = expand(1/2* g'*D*g)
Load data
load data1.mat
b_val = 1000;
xvec = g(:,1);
yvec = g(:,2);
zvec = g(:,3);
Svec = reshape(S, [], 64);
S0vec = reshape(S0, [], 1);
Construct A and b
A = [-xvec.^2/2, -xvec.*yvec, -xvec.*zvec, -yvec.*2/2, -yvec.*zvec, -zvec.^2/2];
b = log(Svec./S0vec);
c = zeros(6, size(S0vec, 1));
for i = 1:size(S0vec, 1)
c(:,i) = A \ b(i,:)';
end
% Tensor for each voxel
D = zeros(3, 3, size(S0vec, 1));
for i = 1:size(S0vec, 1)
D(:,:,i) = [c(1,i), c(2,i), c(3,i); c(2,i), c(4,i), c(5,i); c(3,i), c(5,i), c(6,i)];
end
D;
0 Comments
Answers (0)
See Also
Categories
Find more on Polymers in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!