Replace unwanted voxels with nearest neighbor labeled values
5 次查看(过去 30 天)
显示 更早的评论
banikr
2020-5-3
Hello Matlab experts,
I am trying to remove unwanted electrode voxels from 3D labeled images by replacing the electrodes with nearest skin, bone or other labeled tissue voxels. I experimented with the function here: https://www.mathworks.com/matlabcentral/fileexchange/24723-nearest_neighbour_3d
But that didn't generate the result I expected. The function takes given points(electrodes) and their minimum distant nearest candidate points(other tissue voxels). It considers Euclidean distance.
%sz is the size of the 3D label=> lab
[cxx, cyy, czz] = ind2sub(sz, find(lab==10|lab==4)); % candidate points from either skin(10) or bones(4) or other tissues
[gxx, gyy, gzz] = ind2sub(sz, find(lab==20)); % given poinst from electrodes(20)
[nxx, nyy, nzz] = ind2sub(sz, compute_nearest_neighbour([cxx,cyy,czz], [gxx,gyy,gzz])); % nearest points
lab(gxx,gyy,gzz) = lab(nxx,nyy,nzz);
Any other ways?
Thanks in advance.
采纳的回答
darova
2020-5-15
Solution
lab = lab_superior; % the data shared
elec = zeros(size(lab));
elec(lab==20)=1; % select electrodes
BB = lab;
BB(BB==20) = 0; % remove electrodes
nb = unique(BB); % unique regions
nb(nb==0) = []; % remove '0' from list of regions
[EE,ne] = bwlabeln(elec);
II = BB*0;
for i = 1:ne
E1 = EE == i;
E2 = imdilate(E1,strel('sphere', 4));
for j = nb(:)'
B1 = BB == j; % select one region
tmp = B1 & E2; % compare electrode and region
if any(tmp(:)) % if region and electrode are close
II = II+j*E1;
break
end
end
end
A = II + BB;
Original After processing
更多回答(1 个)
Image Analyst
2020-5-3
The position of x and y are swapped, which makes a difference if the volume is not equal lengths in the x and y direction. Remember ind2sub() returns (row, column) which is NOT (x, y) -- it's (y, x).
36 个评论
banikr
2020-5-3
I see. So the code should be as below?
lab= lab_skin_done;
[czz, cyy, cxx] = ind2sub(sz, find(lab~=20&lab~=0));
[gzz, gyy, gxx] = ind2sub(sz, find(lab==20)); % given poinst from electrodes
[nzz, nyy, nxx] = ind2sub(sz, compute_nearest_neighbour([cxx,cyy,czz], [gxx,gyy,gzz])); % nearest points
lab(gxx,gyy,gzz) = lab(nxx,nyy,nzz);
Is that what you mean?
Input to the function is [xx, yy, zz] sequence right? Opposite order of ind2sub output.
Still I am getting points at unexpected locations in the 3D matrix
banikr
2020-5-5
编辑:banikr
2020-5-5
@darova
You are right. Thanks for the clarification.
sz = size(lab);
[gzz, gxx, gyy] = ind2sub(sz, find(lab==20));
provides the correct voxel indices to the image 'lab'. See the red objects(voxels) in the image below.
Also candidate points are taken accordingly from labeled 4 or 10 voxels.
[czz, cxx, cyy] = ind2sub(sz, find(lab==10|lab==4));
but the function output is not giving me desired output nearest voxels.
[nzz, nxx, nyy] = ind2sub(sz, compute_nearest_neighbour([czz,cxx,cyy], [gzz,gxx, gyy]));
which is nowhere near g(given) values.
>> max([gzz, gxx, gyy])
ans =
129 203 159
>> max([nzz, nxx, nyy])
ans =
170 254 30
Also the maximum values of n values are far from the g ones.
I guess the function is not working...
Any help appreciated.
banikr
2020-5-5
@darova,
How separating each region may help? In fact I wanted to replace the 4 object voxels with nearest neighboring labels.
bwlabeln labels the connected components. I have already the four object voxels labelled as 20.
[gzz, gxx, gyy] = ind2sub(sz, find(lab==20));
Can you explain a bit more?
darova
2020-5-5
Can we start from the start? Can you explain more what you are trying to achieve?
Is this correct?
Are you trying to remove electrodes?
banikr
2020-5-5
Yes, the white region is electrode(lab==20) in a 3D segmented image(170x256x256). I want to replace those electrode voxels with bones(lab==4) or skin(lab==10) or any other tissue labels.
The background is black with lab==0 so I want to avoid that.
banikr
2020-5-5
I wish it was that simple.
some electrodes are embedded inside the skull and there are electrodes objects inside deep brain. So can not wrap up with one label in all the voxels.
darova
2020-5-5
Can you explain more? What electrodes you want to fill/remove and how are you detecting them
banikr
2020-5-5
The 3D label image has 10 different labels. Tissues and electrode.
>> unique(lab)
ans =
0
1
2
3
4
8
10
20
56
20 for electrode and 4 for bones and 10 for skin and others.
The white voxels in the above axial slice are the electrodes.
[gzz, gxx, gyy] = ind2sub(sz, find(lab==20)); % given points from electrodes
dem_lab = zeros(size(lab));
for i=1:numel(gxx)
dem_lab(gzz(i),gxx(i),gyy(i))=1;
end
volumeViewer(dem_lab)
This chunk of code gives the following only electrode voxels.
or simply:
lab(lab~=20)=0;
will provide the same result too.
I want to replace the electrodes with neighboring non-electrode labels. I was thinking about nearest neighbor algorithm. I can hard code it but that will be both memory and time consuming. Trying to find simple or efficient matlab function.
I have segmented the electrodes from CT images which has individual intensity for metals but that is another story.
Please let me know if that is not clear.
Thanks
darova
2020-5-5
You can dilate each electrode and try to merge them with each area/region separately. If there are commond nodes - close region
if any(A1 & A2)
disp('electrode is close to region')
end
See this idea
banikr
2020-5-6
编辑:banikr
2020-5-6
@darova
I think there has been a misunderstanding. The idea you've shared does not GET RID OF the electrode voxels, does it?
Can you clarify how dilation of morphological operations can help? I want to replace the electrodes with nearest labeled voxels not with background or dilated electrode region.
Thanks.
Image Analyst
2020-5-6
编辑:Image Analyst
2020-5-6
I've been sitting out for a few reasons. (1) darova seems to be helping you fine. (2) it looks like it might take more time that I can afford to devote, especially if I have to invent data. And mostly, (3) I have no data to work with. If you had uploaded a small 3-D image in a .mat file that illustrates the need, I might have done more on this. Any way you could upload a small chunk of your volume, say a 500x500x300 labeled subvolume that has all 3 or 4 classes in it?
Image Analyst
2020-5-6
Yes, but it would be easier if you had data. And start your own answer so you can get credit for it if he Accepts or Votes for your Answer. Let's promote you from a "Rising Star" to an MVP.
banikr
2020-5-6
@Image Analyst
Yeah, I usually upload the .mat file but the medical image files have some confidentiality restraints. But I will share small part of it that relates to the problem. Please see the link below:
banikr
2020-5-7
Hi @darova,
Could you explain why did you take the inverted binary?
I1 = ~im2bw(I);
In my case, I is an image with 10 different labels.
lab= lab_superior;
elec = zeros(size(lab));
elec(lab==20)=1;
% [BB,nb] = bwlabeln((lab~=20&lab~=0)&~elec); % label regions?
[BB,nb] = bwlabeln((lab~=0)&~elec); % region without background and without electrodes?
[EE,ne] = bwlabeln(elec); % label all electrode objects.
II = BB*0;
for i = 1:ne
E1 = EE == i; % select one electrode
E2 = imdilate(E1,strel('sphere', 4));
for j = 1:nb
B1 = BB == j; % select one region
% imshowpair(E2,E1|B1) % 3D volumes won't work
% pause(1)
tmp = B1 & E2; % compare electrode and region
if any(tmp(:)) % if region and electrode are close
II = II + j*(B1 | E2); % merge region and electrode
% imshow(II) % for 3D?
% pause(1)
break
end
end
end
In the end, I get II with 3 unique labels. I changed the dilation considering 3D volume. And not clear about this line here:
II = II + j*(B1 | E2);
darova
2020-5-8
- Could you explain why did you take the inverted binary?
- I1 = ~im2bw(I);
It's for bwselect and bwlabel. They want regions to be 'ones' not 'zeros'
- And not clear about this line here:
- II = II + j*(B1 | E2);
You are right, i think better be
II = II + j*E1; % fill not dilated electrode with the same value as region
Did it work? Do you have more questions?
banikr
2020-5-13
Hi @darova,
lab = lab_superior; % the data shared
elec = zeros(size(lab));
elec(lab==20)=1;
I1 = lab;
[BB,nb] = bwlabeln(I1&~elec);
[EE,ne] = bwlabeln(elec);
II = BB*0;
for i = 1:ne
E1 = EE == i;
E2 = imdilate(E1,strel('sphere', 4));
for j = 1:nb
B1 = BB == j; % select one region
tmp = B1 & E2; % compare electrode and region
if any(tmp(:)) % if region and electrode are close
II = II+j*E1;
break
end
end
end
Executing this code gives 'II' with 0 and 1 binary values. The same as elec(with same 475) electrode voxels.
So it didn't work.
Another reason for nearest neighbor is some electrodes share locality with different labels.
So if we dilute the electrodes and assign them common region, will it take the whole electrode as single neighboring label?
I have shared the data in previous comment, were you able to take a look at that?
banikr
2020-5-13
Check this .nii.gz format of the same file:
https://drive.google.com/open?id=1K2GDbycS45DG6ZBw2ecGgGIBq_K2wvIY
darova
2020-5-13
An error occurs
Error using load
Unknown text on line number 1 of ASCII file lab_superior.nii
"\".
What format is it? .nii
banikr
2020-5-13
Yes it is saved in nifti file format.
You can access the file as below:
lab = niftiread('lab_superior.nii.gz');
% volumeViewer(lab) % works
>> unique(lab)
ans =
0
1
2
3
4
8
10
20
darova
2020-5-13
It's because all electrodes merge with BB==1
try
A = BB + II; % eletrodes and brain together
volumeViewer(A)
banikr
2020-5-14
So I include the following:
A = BB + II;
after that, it gives me an image with 28 separate regions. How do I get original label image from that with one less label(20 for electrodes).
>> numel(unique(BB)) == numel(unique(A))
ans =
logical
1
darova
2020-5-14
- Did you run with the code I mentioned here in Previous comment ?
Of course
- How do I get original label image from that with one less label(20 for electrodes).
You didn't tell that you want original labeled image. Just replace BB with your original image (without electrodes)
banikr
2020-5-14
编辑:banikr
2020-5-14
Ahhh, I mentioned couple of times that I want the electrode voxels replaced with neighboring non-electrode labeled voxels.
If we can do that we have original labeled image without electrodes. If there was a misunderstanding sorry for that.
And could you demonstrate with the codes since you have the data too, so far words are very confusing so I post codes which are easy to understand... for example which BB to replace with which lab (do 'without electrodes' in original image means zeroing out lab==20 voxels?) and which line in the code?
darova
2020-5-14
Try this
lab = lab_superior; % the data shared
elec = zeros(size(lab));
elec(lab==20)=1; % select electrodes
BB = lab;
BB(BB==20) = 0; % remove electrodes
nb = unique(BB); % unique regions
[EE,ne] = bwlabeln(elec);
II = BB*0;
for i = 1:ne
E1 = EE == i;
E2 = imdilate(E1,strel('sphere', 4));
for j = nb(:)'
B1 = BB == j; % select one region
tmp = B1 & E2; % compare electrode and region
if any(tmp(:)) % if region and electrode are close
II = II+j*E1;
break
end
end
end
A = II + BB;
banikr
2020-5-14
Hi @darova,
If A ought to be the final result then it removed the electrodes but didn't replace them with neighboring tissues.
lab with electrodes(white voxels)
Thanks for the support, I think I will try some other applications.
darova
2020-5-14
Im sorry, forgot about '0'. The script merged electrodes with background (zero value)
nb = unique(BB); % unique regions
nb(nb==0) = []; % remove '0' from list of regions
Original
After processing
banikr
2020-5-15
@darova,
I think it worked.
The electrode voxels were replaced with mostly skin(brown) voxels.
Thanks for the support!
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Denoising and Compression 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!发生错误
由于页面发生更改,无法完成操作。请重新加载页面以查看其更新后的状态。
您也可以从以下列表中选择网站:
如何获得最佳网站性能
选择中国网站(中文或英文)以获得最佳网站性能。其他 MathWorks 国家/地区网站并未针对您所在位置的访问进行优化。
美洲
- América Latina (Español)
- Canada (English)
- United States (English)
欧洲
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom (English)
亚太
- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)