Extracting Vessel Dimensions from 2D Fluorescence Microscopy Images in MATLAB

11 次查看（过去 30 天）
Alessio
Alessio 2025-3-26
评论： William Rose 2025-3-27
I am working with fluorescence microscopy images of blood vessels in the brain, stored as HDF5 files (.h5). My goal is to extract vessel dimensions (diameters) from each frame in a time series and analyze how they change over time.
Currently, I am using MATLAB to:
  1. Load the .h5 file and extract individual frames.
  2. Apply Gaussian filtering and contrast enhancement to improve vessel visibility.
  3. Perform edge detection (Canny) and skeletonization to extract vessel centerlines.
  4. Compute vessel diameters using the Euclidean distance transform from the skeleton.
However, I am running into several issues:
  • Edge detection is inconsistent and picks up artifacts inside the vessel.
  • The skeleton sometimes extends beyond the vessels or is not continuous.
  • Extracted diameters appear constant, even when the vessel visibly changes size.
Current Approach in MATLAB
Here’s the code I am using to manually select vessels and compute diameters:
vars = dloc.saved_vars;
img = dloc.load_var("reference_img");
img = squeeze(img(:,:,1));
frame_idx = 1000;
frame = double(h5read(dloc.path, '/recording', [1,1,1,frame_idx], [512, 512,1,1]));
frame_smooth = imgaussfilt(frame, 2);
% Create a normalized & contrast-enhanced frame for Canny
frame_enhanced = mat2gray(frame_smooth);
frame_enhanced = imadjust(frame_enhanced, stretchlim(frame_enhanced, [0.01, 0.99]));
%% --- Frangi Vesselness Filter ---
[I_vessel, ~] = FrangiFilter2D(frame_enhanced, struct( ...
'FrangiScaleRange', [1 5], ...
'FrangiScaleRatio', 2, ...
'FrangiBetaOne', 0.5, ...
'FrangiBetaTwo', 15, ...
'BlackWhite', true));
bw_vessel = I_vessel > 0.1;
bw_vessel = bwareaopen(bw_vessel, 100);
bw_vessel = imclose(bw_vessel, strel('disk', 3));
skeleton = bwskel(bw_vessel);
skeleton = bwmorph(skeleton, 'spur', 30);
distance_map = bwdist(~bw_vessel);
diameter_map = distance_map * 2;
% Extract vessel diameters along skeleton
[y, x] = find(skeleton);
diameters = arrayfun(@(i) diameter_map(y(i), x(i)), 1:length(x));
%% --- Edge Detection for Overlay ---
edge_thresholds = [0.1, 0.4];
edges = edge(frame_enhanced, 'Canny', edge_thresholds);
%% --- Visualization: Frangi + Edges + Skeleton + Diameters ---
figure;
imshow(frame_enhanced, []);
hold on;
contour(edges, [1,1], 'g', 'LineWidth', 1);
plot(x, y, 'r.', 'MarkerSize', 5);
scatter(x, y, 15, diameters, 'filled');
colormap(parula);
c = colorbar;
c.Label.String = 'Diameter (pixels)';
title(['Edges + Skeleton + Diameters (Frame ', num2str(frame_idx), ')']);
%% --- Plot Diameter Profile ---
figure;
plot(1:length(diameters), diameters, 'b.-');
xlabel('Skeleton Point Index');
ylabel('Vessel Diameter (pixels)');
title('Vessel Diameter Profile Along Skeleton');
Resulting edges and vessel diameter detection
As you can see, the edges are picked up but not always, and the diameters are calculated on no real basis. They don’t seem to be particularly telling.
My question is:
How can I ensure accurate edge detection that captures the actual vessel boundaries and calculate their diameters?

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

采纳的回答

William Rose
William Rose 2025-3-26
@Alessio,
Consider Better Skeletonization in the File Exchange.
Please upload a frame of the raw .h5 video image.
You write "As you can see, the edges are picked up but not always, and the diameters are calculated on no real basis. They don’t seem to be particularly telling." What you say I can see is not obvious to me. The grayscale and color images in your post do not match either of the figures made by the code you posted. Please post the images that were generated by the code that you posted.
The vessel diameter is a function of both position along the vessel and time: d(x,t), where x is a measure of distance along the axis of the vessel, and t is time (frame number). When you say "Extracted diameters appear constant, even when the vessel visibly changes size", do you mean that d(x) appears constant or d(t) appears constant, or both?
The code you posted does not include a measure of position along the axis of the vessel ("x", in my comment above). This seems like something you may want to do. If you assign an "x" value to each pixel along the centerline, then you will be able to plot diameter versus x, which may be useful.
A number of years ago I wrote code in Labview to extract time-varying vessel diameter from ultrasound videos of human arteries. It worked pretty well. The output, d(t), reqired filtering to remove glitches and gaps, due to occasional changes in insonation (illumination) that occur with ultrasound.
  3 个评论
显示 1更早的评论
William Rose
William Rose 2025-3-27
@Alessio,
If you make a posting that contains mistakes or inconsistencies, please don't edit the original post, if that renders the responses of others non-sensical. We all make mistakes. If you go back and change your post, then people trying to help you cannot figure out what really happened and the discussion makes no sense. The better approach is to post a comment with your revisions. In this case, you went back and changed your code and changed the figures - at least the second figure. In your original posting, the titles and other features of the images did not match the code you posted.
In your original posting, the second image had a skeleton "line" (not a straight line) running down the approximate middle of the central vessel. The skeleton line split at the Y in the vessel, and a branch of the skeleton continued along each of the two vessel branches. Now the skeleton line in the middle is not visible. I think it is because nyour code includes:
plot(x, y, 'r.', 'MarkerSize', 5); % plot skeleton as red dots
scatter(x, y, 15, diameters, 'filled'); % plot circles with diameter 15
The circles plotted in the second line of code above cover up the skeleton. The second line plots circles, all with size 15. The circle colors are given by "diamter" in your code. Is that your intention?
Please post the files and data necessary to allow others to run your code. That includes the file reference_img, and the contents of variable dloc, and the file which is referenced by dloc.path.
William Rose
William Rose 2025-3-27
I see that you extract the red component of the image, with
img = squeeze(img(:,:,1));
That makes sense, if the vessels are labelled with DiI and you use a rhodamine filter, or a similar setup that produces a good signal in the red channel. Otherwise, you could experiment with other approaches such as summing all three colors to make a grayscale image:
img=sum(img,3);
Maybe that will have improved contrast.

请先登录，再进行评论。

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Object Analysis 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by