How to separate into unimodal or bimodal
28 次查看(过去 30 天)
显示 更早的评论
Hello, I have some histograms. How can I separate them through Matlab code and tell which one is unimodal or bimodal ?
thanks in advance
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/934814/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/934819/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/934824/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/934829/image.png)
采纳的回答
Image Analyst
2022-3-22
编辑:Image Analyst
2022-3-22
Maybe use kstest(), lillietest(), or chi2gof().
Or fit the data to two Gaussians and see if the peak locations are separated by enough distance. See attached demos.
Also see links on the right, such as
19 个评论
Jórdan Venâncio Leite
2022-3-23
I'm trying to use the kstest, but i'm not understanding how I'm going to use a standard normal distribution to achieve my goal of differentiating unimodal from bimodal pattern
Jórdan Venâncio Leite
2022-3-24
Yes, however, I still don't know what this part of the code here below does. I know before this is the part where I enter my noisy data. However, after that I got a little confused.
% Convert X and Y into a table, which is the form fitnlm() likes the input data to be in.
tbl = table(X(:), Y(:));
% Define the model as Y = a + b*x + c*exp(-(x-d)^2/e) + d * exp(-(x-f)^2/g)
% Note how this "x" of modelfun is related to big X and big Y.
% x(:, 1) is actually X and x(:, 2) is actually Y - the first and second columns of the table.
modelfun = @(b,x) b(1) + b(2) * x(:, 1) + b(3) * exp(-(x(:, 1) - b(4)).^2/b(5)) + b(6) * exp(-(x(:, 1) - b(7)).^2/b(8));
beta0 = [6, 0.1, 35, 10, 3, 25, 50, 9]; % Guess values to start with. Just make your best guess.
% Now the next line is where the actual model computation is done.
mdl = fitnlm(tbl, modelfun, beta0);
% Now the model creation is done and the coefficients have been determined.
% YAY!!!!
% Extract the coefficient values from the the model object.
% The actual coefficients are in the "Estimate" column of the "Coefficients" table that's part of the mode.
coefficients = mdl.Coefficients{:, 'Estimate'}
% Let's do a fit, but let's get more points on the fit, beyond just the widely spaced training points,
% so that we'll get a much smoother curve.
X = linspace(min(X), max(X), 1920); % Let's use 1920 points, which will fit across an HDTV screen about one sample per pixel.
% Create smoothed/regressed data using the model:
yFitted = coefficients(1) + coefficients(2) * X + coefficients(3) * exp(-(X - coefficients(4)).^2 / coefficients(5)) + ...
coefficients(6) * exp(-(X - coefficients(7)).^2 / coefficients(8));
% Now we're done and we can plot the smooth model as a red line going through the noisy blue markers.
hold on;
plot(X, yFitted, 'r-', 'LineWidth', 2);
grid on;
title('Exponential Regression with fitnlm()', 'FontSize', fontSize);
xlabel('X', 'FontSize', fontSize);
ylabel('Y', 'FontSize', fontSize);
legendHandle = legend('Noisy Y', 'Fitted Y', 'Location', 'northeast');
legendHandle.FontSize = 25;
Jórdan Venâncio Leite
2022-3-24
编辑:Jórdan Venâncio Leite
2022-3-24
My code is like this
load('image66.mat');
img=PerfilCinzasVertical;
N = histogram(img,12);
using the archive attached
Jórdan Venâncio Leite
2022-3-28
Image Analyst, i can keep using this model here in my case which is a histogram? MODEL: Y = a + b*x + c*exp(-(x-d)^2/e) + d * exp(-(x-f)^2/g). How do I convert the histogram into a table to be used in the line 10 of my code knowing that the variable of this histogram in the workspace is a graph and not values?
clc
clear
close all
load('image66.mat');
img=PerfilCinzasVertical;
N = histogram(img,12);
% Convert the histogram into a table
tbl = table(X(:), Y(:));
% Define the model as Y = a + b*x + c*exp(-(x-d)^2/e) + d * exp(-(x-f)^2/g)
% Note how this "x" of modelfun is related to big X and big Y.
% x(:, 1) is actually X and x(:, 2) is actually Y - the first and second columns of the table.
modelfun = @(b,x) b(1) + b(2) * x(:, 1) + b(3) * exp(-(x(:, 1) - b(4)).^2/b(5)) + b(6) * exp(-(x(:, 1) - b(7)).^2/b(8));
beta0 = [6, 0.1, 35, 10, 3, 25, 50, 9]; % Guess values to start with. Just make your best guess.
% Now the next line is where the actual model computation is done.
mdl = fitnlm(tbl, modelfun, beta0);
% Now the model creation is done and the coefficients have been determined.
% YAY!!!!
% Extract the coefficient values from the the model object.
% The actual coefficients are in the "Estimate" column of the "Coefficients" table that's part of the mode.
coefficients = mdl.Coefficients{:, 'Estimate'}
% Let's do a fit, but let's get more points on the fit, beyond just the widely spaced training points,
% so that we'll get a much smoother curve.
X = linspace(min(X), max(X), 1920); % Let's use 1920 points, which will fit across an HDTV screen about one sample per pixel.
% Create smoothed/regressed data using the model:
yFitted = coefficients(1) + coefficients(2) * X + coefficients(3) * exp(-(X - coefficients(4)).^2 / coefficients(5)) + ...
coefficients(6) * exp(-(X - coefficients(7)).^2 / coefficients(8));
% Now we're done and we can plot the smooth model as a red line going through the noisy blue markers.
hold on;
plot(X, yFitted, 'r-', 'LineWidth', 2);
grid on;
title('Exponential Regression with fitnlm()', 'FontSize', fontSize);
xlabel('X', 'FontSize', fontSize);
ylabel('Y', 'FontSize', fontSize);
legendHandle = legend('Noisy Y', 'Fitted Y', 'Location', 'northeast');
legendHandle.FontSize = 25;
% Set up figure properties:
% Enlarge figure to full screen.
set(gcf, 'Units', 'Normalized', 'OuterPosition', [0 0 1 1]);
% Get rid of tool bar and pulldown menus that are along top of figure.
% set(gcf, 'Toolbar', 'none', 'Menu', 'none');
% Give a name to the title bar.
set(gcf, 'Name', 'Demo by ImageAnalyst', 'NumberTitle', 'Off')
Jórdan Venâncio Leite
2022-3-29
编辑:Jórdan Venâncio Leite
2022-3-29
Image Analyst, which model you suggest me to use with the figure 2 (second histogram)?
Thanks in advance
Image Analyst
2022-3-30
Looks more of a log normal or exponential decay for the left hump and then possibly a gaussian for the right hump.
Jórdan Venâncio Leite
2022-4-8
Image Analyst. I'm having doubts about which model to use to represent the gaussian. Can you give me some tips ?
clc;
clear
close all
load('image66.mat');
img=PerfilCinzasVertical;
N = histogram(img,12);
% Get the bin values into their own separate variable.
x = N.Values
% Get the bin counts into their own separate variable.
y = N.BinEdges(1:12)
% Convert the histogram into a table
tbl = table(x(:), y(:));
% Now the next line is where the actual model computation is done.
mdl = fitnlm(tbl, modelfun, beta0);
% Now the model creation is done and the coefficients have been determined.
% YAY!!!!
% Extract the coefficient values from the the model object.
% The actual coefficients are in the "Estimate" column of the "Coefficients" table that's part of the mode.
coefficients = mdl.Coefficients{:, 'Estimate'}
% Let's do a fit, but let's get more points on the fit, beyond just the widely spaced training points,
% so that we'll get a much smoother curve.
X = linspace(min(X), max(X), 1920); % Let's use 1920 points, which will fit across an HDTV screen about one sample per pixel.
% Create smoothed/regressed data using the model:
yFitted = coefficients(1) + coefficients(2) * X + coefficients(3) * exp(-(X - coefficients(4)).^2 / coefficients(5)) + ...
coefficients(6) * exp(-(X - coefficients(7)).^2 / coefficients(8));
% Now we're done and we can plot the smooth model as a red line going through the noisy blue markers.
hold on;
plot(X, yFitted, 'r-', 'LineWidth', 2);
grid on;
title('Exponential Regression with fitnlm()', 'FontSize', fontSize);
xlabel('X', 'FontSize', fontSize);
ylabel('Y', 'FontSize', fontSize);
legendHandle = legend('Noisy Y', 'Fitted Y', 'Location', 'northeast');
legendHandle.FontSize = 25;
% Set up figure properties:
% Enlarge figure to full screen.
set(gcf, 'Units', 'Normalized', 'OuterPosition', [0 0 1 1]);
% Get rid of tool bar and pulldown menus that are along top of figure.
% set(gcf, 'Toolbar', 'none', 'Menu', 'none');
% Give a name to the title bar.
set(gcf, 'Name', 'Demo by ImageAnalyst', 'NumberTitle', 'Off')
Image Analyst
2022-4-9
Try this:
clc; % Clear the command window.
close all; % Close all figures (except those of imtool.)
clear; % Erase all existing variables. Or clearvars if you want.
workspace; % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 20;
s = load('image66.mat')
subplot(2, 1, 1);
img = s.PerfilCinzasVertical;
plot(img, 'b-');
grid on;
title('Original Data', 'FontSize', fontSize);
subplot(2, 1, 2);
histogramObject = histogram(img, 12)
xlabel('Value', 'FontSize', fontSize);
ylabel('Count', 'FontSize', fontSize);
grid on;
% Convert the histogram into a table
X = histogramObject.BinEdges(1:end-1);
Y = histogramObject.BinCounts;
tbl = table(X(:), Y(:))
% Define the model as Y = a + b*x + c*exp(-(x-d)^2/e) + d * exp(-(x-f)^2/g)
% Note how this "x" of modelfun is related to big X and big Y.
% x(:, 1) is actually X and x(:, 2) is actually Y - the first and second columns of the table.
modelfun = @(b,x) b(1) * exp(-(x(:, 1) - b(2)).^2/b(3)) + b(4) * exp(-(x(:, 1) - b(5)).^2/b(6));
beta0 = [100, 190, 200, 150, 235, 200]; % Guess values to start with. Just make your best guess.
% Plot guess
xl = xlim; % Get range of x axis.
xGuess = linspace(xl(1), xl(2), 500)';
yGuess = modelfun(beta0, xGuess)
hold on;
plot(xGuess, yGuess, 'r-');
hold off;
% Now the next line is where the actual model computation is done.
mdl = fitnlm(tbl, modelfun, beta0);
% Now the model creation is done and the coefficients have been determined.
% YAY!!!!
% Extract the coefficient values from the the model object.
% The actual coefficients are in the "Estimate" column of the "Coefficients" table that's part of the mode.
coefficients = mdl.Coefficients{:, 'Estimate'}
% Let's do a fit, but let's get more points on the fit, beyond just the widely spaced training points,
% so that we'll get a much smoother curve.
X = linspace(min(X), max(X), 1920); % Let's use 1920 points, which will fit across an HDTV screen about one sample per pixel.
% Create smoothed/regressed data using the model:
yFitted = coefficients(1) * exp(-(X - coefficients(2)).^2 / coefficients(3)) + ...
coefficients(4) * exp(-(X - coefficients(5)).^2 / coefficients(6));
% Now we're done and we can plot the smooth model as a red line going through the noisy blue markers.
hold on;
plot(X, yFitted, 'r-', 'LineWidth', 2);
grid on;
title('Exponential Regression with fitnlm()', 'FontSize', fontSize);
xlabel('X', 'FontSize', fontSize);
ylabel('Y', 'FontSize', fontSize);
legendHandle = legend('Noisy Y', 'Guess', 'Fitted Y', 'Location', 'northwest');
legendHandle.FontSize = 25;
% Set up figure properties:
% Enlarge figure to full screen.
set(gcf, 'Units', 'Normalized', 'OuterPosition', [0 0 1 1]);
% Get rid of tool bar and pulldown menus that are along top of figure.
% set(gcf, 'Toolbar', 'none', 'Menu', 'none');
% Give a name to the title bar.
set(gcf, 'Name', 'Demo by ImageAnalyst', 'NumberTitle', 'Off')
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/958075/image.png)
Jórdan Venâncio Leite
2022-4-9
Thank You Very Much for your answer ImageAnalyst !!! Do you think I can have automatic processing of my histograms? That is, for each different histogram, do you believe I can automatically determine if it is unimodal, bimodal or multimodal? or do you believe that i would need to manually configure a model every time i have a different histogram? thanks in advance!
Image Analyst
2022-4-9
Yes, I'd probably first fit a single mode gaussian, and then a two Gaussian model. Then whichever one has the lowest MSE I'd say is the one it is.
Image Analyst
2022-4-9
Mean Square Error. You can use immse() or just use the formula
mse = mean((trueValue - estimatedValue) .^ 2)
Jórdan Venâncio Leite
2022-4-9
Image Analyst, are you saying that it is not possible to make this distinction in the histograms with respect to whether it is unimodal, bimodal or multimodal automatically? txs in advance.
Image Analyst
2022-4-9
To be clear, I am not a statistician. I don't know that much about kstest, lillietest, anova, and other stats things. I suggest you ask your local statistician.
Why do you need to know anyway? What would you do with your signal if you knew the histogram indicated your signal came from two populations or one population?
Jórdan Venâncio Leite
2022-4-9
I just need to know if the histogram is unimodal, bimodal or multimodal. This is my only problem for now.
Image Analyst
2022-4-9
Strange. So you don't know why you need to know? Is it your homework then, and it's just something the professor asked? Did you discuss some method in class? Or if it's not your homework, why don't we just forget about that and figure out what to do with your original signal. Of course any signal that has even a little bit of noise could be considered the composition of two Gaussians. So let's just assume it is. OK, now what? What do you now need to do with your original signal knowing that its histogram is two Gaussians?
Well sorry I couldn't help more (I'm not an expert statistician) but good luck.
Jórdan Venâncio Leite
2022-4-10
THANKS for your answer IMAGE ANALYST!! See if you can help me with this https://www.mathworks.com/matlabcentral/answers/1693350-how-to-define-if-a-histogram-is-unimodal-or-bimodal. It's almost finished already .
更多回答(0 个)
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Probability Distributions 的更多信息
标签
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 (한국어)