Representing three variables in a three dimension plot

Question

I have a problem dealing with 3rd dimension plot for three variables.

I have three matrices: Temperature, Humidity and Power. During one year, at every hour, each one of the above were measured. So, we have for each matrix 365*24 = 8760 points. Then, one average point is taken every day. So,

Tavg = 365 X 1
Havg = 365 X 1
Pavg = 365 X 1

In electrical point of veiw, the power depends on the temperature and humidity. I want to discover this relation using a three dimensional plot.

I tried using mesh, meshz, surf, plot3, and many other commands in MATLAB but unfortunately I couldn't get what I want. For example, let us take first 10 days. Here, every day is represented by average temperature, average humidity and average power.

Tavg = [18.6275
   17.7386
   15.4330
   15.4404
   16.4487
   17.4735
   19.4582
   20.6670
   19.8246
   16.4810];

Havg = [75.7105
   65.0892
   40.7025
   45.5119
   47.9225
   62.8814
   48.1127
   62.1248
   73.0119
   60.4168];

Pavg = [13.0921
   13.7083
   13.4703
   13.7500
   13.7023
   10.6311
   13.5000
   12.6250
   13.7083
   12.9286];

How do I represent these matrices by three dimension plot?

Answer 1

The challenge is that the 3-D surface plotting functions (mesh, surf, etc.) are looking for a 2-D matrix of z values. So to use them you need to construct such a matrix from the data.

Currently the data is sea of points in 3-D space, so, you have to map these points to a surface. A simple approach to this is to divide up the X-Y (temperature-humidity) plane into bins and then take the average of all of the Z (power) data. Here is some sample code for this that uses accumarray() to compute the averages for each bin: % Specify bin sizes Tbin = 3; Hbin = 20;

% Create binned average array
% First create a two column array of bin indexes to use as subscripts
subs = [round(Havg/Hbin)+1, round(Tavg/Tbin)+1];

% Now create the Z (power) estimate as the average value in each bin
Pest = accumarray(subs,Pavg,[],@mean);

% And the corresponding X (temp) & Y (humidity) vectors
Tval = Tbin/2:Tbin:size(Pest,2)*Tbin;
Hval = Hbin/2:Hbin:size(Pest,1)*Hbin;

% And create the plot
figure(1)
surf(Tval, Hval, Pest)
xlabel('Temperature')
ylabel('Humidity')
zlabel('Power')
title('Simple binned average')

xlim([14 24])
ylim([40 80])

The graph is a bit coarse (can't post image yet, since I am new) because we only have a few data points. We can enhance the visualization by removing any empty bins by setting their value to NaN. Also the binning approach hides any variation in the Z (power) data so we can also overlay the orgional point cloud using plot3 without drawing connecting lines. (Again no image b/c I am new)

Additional code for the final plot: %% Expanded Plot

% Remove zeros (useful with enough valid data)
%Pest(Pest == 0) = NaN;

% First the original points
figure(2)
plot3(Tavg, Havg, Pavg, '.')

hold on
% And now our estimate
% The use of 'FaceColor' 'Interp' uses colors that "bleed" down the face
% rather than only coloring the faces away from the origin
surfc(Tval, Hval, Pest, 'FaceColor', 'Interp')

% Make this plot semi-transparent to see the original dots anb back side
alpha(0.5)

xlabel('Temperature')
ylabel('Humidity')
zlabel('Power')
grid on

title('Nicer binned average')
xlim([14 24])
ylim([40 80])

Answer 2

I think you're asking for a surface fit for your data. The Curve Fitting Toolbox handles this nicely:

% Fit model to data.
ft = fittype( 'poly11' );
fitresult = fit( [Tavg, Havg], Pavg, ft);

% Plot fit with data.
plot( fitresult, [xData, yData], zData );
legend( 'fit 1', 'Pavg vs. Tavg, Havg', 'Location', 'NorthEast' );
xlabel( 'Tavg' );
ylabel( 'Havg' );
zlabel( 'Pavg' );
grid on

If you don't have the Curve Fitting Toolbox, you can use the backslash operator:

% Find the coefficients.
const = ones(size(Tavg));
coeff = [Tavg Havg const] \ Pavg; 

% Plot the original data points
clf
plot3(Tavg,Havg,Pavg,'r.','MarkerSize',20);
hold on

% Plot the surface.
[xx, yy] = meshgrid( ...
    linspace(min(Tavg),max(Tavg)) , ...
    linspace(min(Havg),max(Havg)) ); 
zz = coeff(1) * xx + coeff(2) * yy + coeff(3);
surf(xx,yy,zz)
title(sprintf('z=(%f)*x+(%f)*y+(%f)',coeff))
grid on
axis tight

Both of these fit a linear polynomial surface, i.e. a plane, but you'll probably want to use something more complicated. Both of these techniques can be adapted to this situation. There's more information on this subject at mathworks.com: How can I determine the equation of the best-fit line, plane, or N-D surface using MATLAB?.

Answer 3

You might want to look at Delaunay triangulation:

tri = delaunay(Tavg, Havg);
trisurf(tri, Tavg, Havg, Pavg);

Using your example data, this code generates an interesting 'surface'. But I believe this is another way of doing what you want.

Answer 4

You might also try the GridFit tool by John D'Errico from MATLAB Central. This tool produces a surface similar to interpolating between the data points (as is done by MATLAB's griddata) but with cleaner results because it smooths the resulting surface. Conceptually multiple datapoints for nearby or overlapping X,Y coordinates are averaged to produce a smooth result rather than noisy "ripples." The tool also allows for some extrapolation beyond the data points. Here is a code example (assuming the GridFit Tool has already been installed):

%Establish points for surface
num_points = 20;
Tval = linspace(min(Tavg),max(Tavg),num_points);
Hval = linspace(min(Havg),max(Havg),num_points);

%Do the fancy fitting with smoothing
Pest = gridfit(Tavg, Havg, Pavg, Tval, Hval);

%Plot results
figure(5)
surfc(XI,YI,Pest, 'FaceColor', 'Interp')

To produce an even nicer plot, you can add labels, some transparancy and overlay the original points:

alpha(0.5)
hold on
plot3(Tavg,Havg,Pavg,'.')

xlabel('Temperature')
ylabel('Humidity')
zlabel('Power')
grid on

title('GridFit')

PS: @upperBound: Thanks for the Delaunay triangulation tip. That seems like the way to go if you want to go through each of the points. I am a newbie so can't comment yet.

Answer 5

Below is your solution:

1. Save/write the Myplot3D function

   function [x,y,V]=Myplot3D(X,Y,Z)
   x=linspace(X(1),X(end),100);
   y=linspace(Y(1),Y(end),100);
   [Xt,Yt]=meshgrid(x,y);
   V=griddata(X,Y,Z,Xt,Yt);
   

2. Call the following from your command line (or script)

   [Tavg_new,Pavg_new,V]=Myplot3D(Tavg,Pavg,Havg);
   surf(Tavg_new,Pavg_new,V)
   colormap jet;
   xlabel('Temperature')
   ylabel('Power/Pressure')
   zlabel('Humidity')