You are now following this question
- You will see updates in your followed content feed.
- You may receive emails, depending on your communication preferences.
How to separate into unimodal or bimodal
28 views (last 30 days)
Show older comments
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)
Accepted Answer
Image Analyst
on 22 Mar 2022
Edited: Image Analyst
on 22 Mar 2022
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 Comments
Jórdan Venâncio Leite
on 23 Mar 2022
Thanks for the answer Image Analyst!
Jórdan Venâncio Leite
on 23 Mar 2022
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
Image Analyst
on 23 Mar 2022
Did you try fitting the data to two Gaussians using the other demo?
Jórdan Venâncio Leite
on 24 Mar 2022
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
on 24 Mar 2022
Edited: Jórdan Venâncio Leite
on 24 Mar 2022
My code is like this
load('image66.mat');
img=PerfilCinzasVertical;
N = histogram(img,12);
using the archive attached
Jórdan Venâncio Leite
on 28 Mar 2022
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')
Image Analyst
on 28 Mar 2022
X would be the bin values, and Y would be the bin counts.
Jórdan Venâncio Leite
on 29 Mar 2022
Edited: Jórdan Venâncio Leite
on 29 Mar 2022
Image Analyst, which model you suggest me to use with the figure 2 (second histogram)?
Thanks in advance
Image Analyst
on 30 Mar 2022
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
on 8 Apr 2022
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
on 9 Apr 2022
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
on 9 Apr 2022
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
on 9 Apr 2022
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
on 9 Apr 2022
Mean Square Error. You can use immse() or just use the formula
mse = mean((trueValue - estimatedValue) .^ 2)
Jórdan Venâncio Leite
on 9 Apr 2022
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
on 9 Apr 2022
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
on 9 Apr 2022
I just need to know if the histogram is unimodal, bimodal or multimodal. This is my only problem for now.
Image Analyst
on 9 Apr 2022
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
on 10 Apr 2022
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 .
More Answers (0)
See Also
Categories
Find more on Probability Distributions in Help Center and File Exchange
Tags
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!An Error Occurred
Unable to complete the action because of changes made to the page. Reload the page to see its updated state.
Select a Web Site
Choose a web site to get translated content where available and see local events and offers. Based on your location, we recommend that you select: .
You can also select a web site from the following list
How to Get Best Site Performance
Select the China site (in Chinese or English) for best site performance. Other MathWorks country sites are not optimized for visits from your location.
Americas
- América Latina (Español)
- Canada (English)
- United States (English)
Europe
- 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)
Asia Pacific
- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)