Bivariate histogram without hist3

Hello, as illustrated in the image below I'm trying to create an bivariate histogram plot. But I can't use the hist3 function (don't have the statistic toolbox).

Is there a way to do this with Matlab 2012b without the hist3 function?

Answer 1

[Adam Danz]

There are multiple ways. histrogram2() would be the best alternative but it wasn't introduced until r2015b. I recreated the image you attached by using heatmap() but it could also be done using imagesc() or other methods. Look out for some small difference such as how the bins are labeled.

After loading the data I remove NaNs and then bin the data using histcounts(). Unfortunately heatmap() cannot label the edges of the bins; it only labels their centers. So I assign the data to the bin centers. The heatmap() function does the rest.

load carbig

% Remove NaNs
nanIdx = isnan(MPG) | isnan(Weight); 
MPG(nanIdx) = [];
Weight(nanIdx) = []; 

% bin the data
nBins = 10; %number of bins (there will be nBins + 1 edges)
[~, mpgEdge, mpgBin] = histcounts(MPG, nBins);
[~, wgtEdge, wgtBin] = histcounts(Weight, nBins);

% Calculate center of each bin
mpgBinCnt = mpgEdge(2:end) - (diff(mpgEdge)/2); 
wgtBinCnt = wgtEdge(2:end) - (diff(wgtEdge)/2); 

% Assign data to bins
MPGb = mpgBinCnt(mpgBin)'; 
Wightb = wgtBinCnt(wgtBin)'; 

% Put into table and plot heatmap
tbl = table(MPGb, Wightb); 
hh = heatmap(tbl, 'MPGb', 'Wightb', 'CellLabelColor', 'none');

% Revers y axis
hh.YDisplayData = flipud(hh.YDisplayData); 

% Label the axes
hh.XLabel = 'MPG'; 
hh.YLabel = 'Weight'; 
hh.Title = 'hist3 simulation'; 

% Change colormap
colormap parula

[leonidas86]

Thanks for your response Adam. Because I'm using Matlab 2012b I can't use the functions histcounts and heatmap.

I will try it with histc and imagesc. It's problematic to get the same length for MPGb and Weightb with histc because there is no input for the bin number..

[leonidas86]

So I tried it with the commands from 2012b as follows:

% bin the data
mpgMax=max(MPG)+1;
[~,mpgID] = histc(MPG, 0:mpgMax);
wgtMax=max(Weight)+1;
[~,wgtID] = histc(Weight, 0:wgtMax);

mpgBinsMax = accumarray(mpgID,MPG,[],@max);
mpgBinsMin = accumarray(mpgID,MPG,[],@min);
mpgDiff = mpgBinsMax - mpgBinsMin;
mpgCenter = mpgBinsMin + mpgDiff/2;

wgtBinsMax = accumarray(wgtID,MPG,[],@max);
wgtBinsMin = accumarray(wgtID,MPG,[],@min);
wgtDiff = wgtBinsMax - wgtBinsMin;
wgtCenter = wgtBinsMin + wgtDiff/2;

% assign data to bins
MPGb = mpgCenter(mpgID)';
Weightb = wgtCenter(wgtID)';

% Put into table and plot heatmap
histo = horzcat(MPGb, Weightb); 
hh = imagesc(histo);

% Change colormap
colormap jet

But I only get this crazy result:

[Adam Danz]

Oops, I forgot to check which functions weren't available in 2012b.

Here's the updated code using hisc() and imagesc(). imagesc() is more flexible than heatmap() so this version is better, anyway. Note that I haven't checked if any of the other functions I use here were in 2012b but I think we should be good.

imagesc() doesn't frame each block with black lines but you could easily add that in if you wish.

load carbig

% Remove NaNs
nanIdx = isnan(MPG) | isnan(Weight); 
MPG(nanIdx) = [];
Weight(nanIdx) = []; 

% bin the data
nBins = 10; %number of bins (there will be nBins + 1 edges)
mpgEdge = linspace(min(MPG),max(MPG),nBins); 
wgtEdge = linspace(min(Weight),max(Weight),nBins); 
[mpgN, mpgBin] = histc(MPG, mpgEdge);
[wgtN, wgtBin] = histc(Weight, wgtEdge);

% count number of elements per (x,y) pair
[xIdx, yIdx] = meshgrid(1:nBins, 1:nBins); 
xyPairs = [xIdx(:), yIdx(:)]; 

nPerBin = zeros(size(xyPairs,1),1); 
for i = 1:size(xyPairs,1)
    nPerBin(i) = sum(ismember([mpgBin, wgtBin], xyPairs(i,:), 'rows')); 
end

% Reshape nPerBin to grid size
nPerBinGrid = reshape(nPerBin, [nBins, nBins]);

% plot data 
figure
ih = imagesc(nPerBinGrid, [min(nPerBin), max(nPerBin)]);

% Reverse y axis
set(gca, 'YDir', 'normal');

% Change colormap
colormap parula

% Label the axes
xlabel('MPG')
ylabel('Weight')
title('hist3 simulation using matlab 2012');

% set tick marks; label edges
mpgEdge(end+1) = mpgEdge(end) + (mpgEdge(2)-mpgEdge(1));  %add 1 more tick
wgtEdge(end+1) = wgtEdge(end) + (wgtEdge(2)-wgtEdge(1));  %add 1 more tick
set(gca, 'XTick', 0.5: 1 : nBins+0.5, 'XTickLabel', round(mpgEdge*10)/10)
set(gca, 'YTick', 0.5: 1 : nBins+0.5, 'YTickLabel', round(wgtEdge*10)/10)

% add colorbar
colorbar()

[leonidas86]

Thanks. It's working perfect!

[leonidas86]

Hello Adam,

I'm trying to add a plot to my image. I added the lines: ih = imagesc(countPerBinGrid, [min(countPerBin), max(countPerBin)]);

hold on
plot(TestArray(:,1),TestArray(:,2),'-','LineWidth',2,'Color',Red);

But the plot doesn't appear. Where is my fault?

[Adam Danz]

Let me guess, the values of TestArray are between ~9 and ~50 along the x axis and ~1613 to ~5500 along the y axis? These aren't the actual values of the x and y coordinates in the plot.

To understand this, run my example code up to this line but no further.

...
ih = imagesc(nPerBinGrid, [min(nPerBin), max(nPerBin)]);

You'll see that the (x,y) coordinates span from 1:10. These are the actual values along the x,y axes and notice that they are centered in each bin.

Later we flip the Y axis and then assign different x and y tick labels using

set(gca, 'XTick', ..., 'XTickLabel', ...)
set(gca, 'YTick', ..., 'YTickLabel', ...)

Even though the end result spans from 9 to 50.8 along the x axis and 1613 to 5532 along the y axis, those are just labels. The actual coordinates are just the bin indices 1:10. That's the danger of assigning alternative tick labels.

To confirm this, set your TestArray equal to 1:10 and plot the diagonal line across the grid:

TestArray = [(1:10)', (1:10)']; 
hold on
ph = plot(TestArray(:,1),TestArray(:,2), 'k-')

I'll write another comment in a bit with ways to get around this problem.

[Adam Danz]

I adjusted my code so that you no longer have to create XTickLabels and YTickLabels which isn't good practice (my bad!). This version specifies the x and y coordinates in imagesc() and it removes the tick label lines. The only 2 changes were

1. the line that calls imagesc()
2. I removed the set(gca, 'X/YTick', ..., 'X/YTickLabel', ...) lines

This version is (more) correct.

load carbig

% Remove NaNs
nanIdx = isnan(MPG) | isnan(Weight); 
MPG(nanIdx) = [];
Weight(nanIdx) = []; 

% bin the data
nBins = 10; %number of bins (there will be nBins + 1 edges)
mpgEdge = linspace(min(MPG),max(MPG),nBins); 
wgtEdge = linspace(min(Weight),max(Weight),nBins); 
[~, mpgBin] = histc(MPG, mpgEdge);
[~, wgtBin] = histc(Weight, wgtEdge);

% count number of elements per (x,y) pair
[xIdx, yIdx] = meshgrid(1:nBins, 1:nBins); 
xyPairs = [xIdx(:), yIdx(:)]; 

nPerBin = zeros(size(xyPairs,1),1); 
for i = 1:size(xyPairs,1)
    nPerBin(i) = sum(ismember([mpgBin, wgtBin], xyPairs(i,:), 'rows')); 
end

% Reshape nPerBin to grid size
nPerBinGrid = reshape(nPerBin, [nBins, nBins]);

% plot data 
figure
ih = imagesc(mpgEdge, wgtEdge, nPerBinGrid, [min(nPerBin), max(nPerBin)]);  

% Reverse y axis
set(gca, 'YDir', 'normal');

% Change colormap
colormap parula

% Label the axes
xlabel('MPG')
ylabel('Weight')
title('hist3 simulation using matlab 2012');                         

% add colorbar
colorbar()

Now your x,y tick marks will be cleaner but they won't necessarily be centered on each bin (which is generally OK). If you specifically want to label the edges, you can add those lines back in. You can now add additional lines using the actual x and y coordinates that are labeled...

TestArray = [(10:5:45)', (1500:500:5000)'];
hold on
ph = plot(TestArray(:,1),TestArray(:,2), 'k-')

[leonidas86]

Thanks, it works perfect.

One last question: How can I scale the colorbar axis logarithmically?

[Adam Danz]

This will help,

https://www.mathworks.com/matlabcentral/answers/393174-how-can-i-make-a-logarithmic-color-plot-and-colorbar

[leonidas86]

This solution only works with newer Matlab versions but not with Matlab 2012b.

[Adam Danz]

These seem worth a try

use caxis: https://www.mathworks.com/matlabcentral/answers/203697-logarithmic-scale-for-colorbar

older: https://www.mathworks.com/matlabcentral/answers/100066-how-do-i-create-a-logarithmic-scale-colormap-or-colorbar