How to calculate some statistical moments for every cell of a cell array?

1 次查看(过去 30 天)
Hello, I produced pairwise distances from some latitude and longitudes, then I want to calculate normal pdf and logarithmic pdf but I get an error message. My aim is to find if these distances are normally distributed or not. (Distances are produced according to the grid i defined onto a map). You can find the data, the code and the error message below:
clear all;
addpath('D:\Desktop\BOUN\JOB\Slip_Rates\Slip_Data\MAP');
LAT1=39; LAT2=42; LON1=29.0; LON2=41.0; sf = 1;
m_proj('albers equal-area','lat',[LAT1 LAT2],'long',[LON1 LON2],'rect','on');m_gshhs_h('color',[.5 .5 .5]);hold on;
m_grid('linewi',1,'linest','none','tickdir','in','fontsize',10);
all = load('all_velocities.txt'); lon1=all(:,1); lat1=all(:,2); ve1=all(:,3); vn1=all(:,4);
% GRIDDING
dLO = .5*2; dLA = .3636*2;
lon = [LON1:dLO:LON2];
lat = [LAT1:dLA:LAT2];
nlat = length(lat); nlon = length(lon); DIST = cell(nlat, nlon);
for i = 1:nlat
for j = 1:nlon
ind = find(abs(lat1-lat(i))<dLA/2 & (abs(lon1-lon(j))<dLO/2));
DENSITY(i,j) = length(ind); % DENSITY OF THE GRIDS
points = [reshape(lat1(ind),[],1), reshape(lon1(ind),[],1)];
P{i,j} = points; % ELEMENTS IN EACH GRID
DIST{i,j} = pdist(points, 'euclidean'); % PAIRWISE DISTANCE BETWEEN ELEMENTS IN EACH GRID
end
end
for i = 1:nlat
for j = 1:nlon
ind = DIST{i,j};
SDIST(i,j) = std(DIST{i,j}); % STANDARD DEVIATION OF DISTANCES
MDIST(i,j) = mean(DIST{i,j}); % MEAN OF DISTANCES
MEDIST(i,j) = median(DIST{i,j}); % MEDIAN OF DISTANCES
SDIST(isnan(SDIST))=0; MDIST(isnan(MDIST))=0; MEDIST(isnan(MEDIST))=0;
LD(i,j) = exp(MDIST(i,j)+0.5*(SDIST(i,j).^2)); % CALCULATING LOGARITMIC PDF
GD(i,j) = (1/((std(DIST{i,j})*sqrt(2*pi))))*exp(-((ind-mean(DIST{i,j})).^2)/(2*(std(DIST{i,j}).^2))); % CALCULATING GAUSSIAN PDF
end
end
C = abs(MDIST - MEDIST); m_pcolor(lon-dLO/2, lat-dLA/2, C);
m_quiver(lon1,lat1,sf.*ve1,sf.*vn1,1,'w','filled','AutoScale','off','linewidth',1.5); hold off;
Unable to perform assignment because the size of the left side is 1-by-1 and the size of the right side is 1-by-0.
Error in all_velocities4 (line 48)
GD(i,j) = (1/((std(DIST{i,j})*sqrt(2*pi))))*exp(-((ind-mean(DIST{i,j})).^2)/(2*(std(DIST{i,j}).^2))); % CALCULATING GAUSSIAN
PDF
  6 个评论
Guillaume
Guillaume 2020-2-4
In your 2nd set of loops you have:
ind = DIST{i,j};
%...
GD(i,j) = (1/((std(DIST{i,j})*sqrt(2*pi))))*exp(-((ind-mean(DIST{i,j})).^2)/(2*(std(DIST{i,j}).^2))); % CALCULATING GAUSSIAN PDF
%ind is used here ----------------------------------^
What is the point of ind? First, why such a poor variable name which implies indices when it actually contains distances. Secondly, it's exactly the same as dist{i,j} and you use dist{i,j} everywhere but this one instance in the GD calculation. Why? As far as I can tell it just makes it harder to understand the expression and you could just get rid of the variable. i.e have DIST{i,j} = mean(DIST{i,j}) instead of ind - DIST{i, j}.
Unless ind is meant to be something else...

请先登录,再进行评论。

回答(2 个)

Jakob B. Nielsen
Jakob B. Nielsen 2020-2-4
A good way to troubleshoot is to insert a keyboard line - it pauses the code, and allows you to see "what matlab sees" mid-execution.
keyboard
GD(i,j) = (1/((std(DIST{i,j})*sqrt(2*pi))))*exp(-((ind-mean(DIST{i,j})).^2)/(2*(std(DIST{i,j}).^2))); % CALCULATING GAUSSIAN PDF
If you pause it here, when i=1 and j=1 std(DIST{i,j}) gives 0. That means you try to divide by zero, never a good idea. However, matlab prints it as NaN and that goes into your GD variable.
However, when i=1 and j=2, DIST{i,j} is an empty double. So you try to place an empty double (1x0 dimension) into GD(1,2) (a 1,1 position within GD).
The root of your problem seems to be with this part of your code. For i=1, j=2 your ind findes no indexes that live up to your limits. That gives your empty cell which goes on to give your error later.
for i = 1:nlat
for j = 1:nlon
ind = find(abs(lat1-lat(i))<dLA/2 & (abs(lon1-lon(j))<dLO/2));
DENSITY(i,j) = length(ind); % DENSITY OF THE GRIDS
points = [reshape(lat1(ind),[],1), reshape(lon1(ind),[],1)];
P{i,j} = points; % ELEMENTS IN EACH GRID
DIST{i,j} = pdist(points, 'euclidean'); % PAIRWISE DISTANCE BETWEEN ELEMENTS IN EACH GRID
keyboard %I inserted a keyboard here, as well
end
end
  5 个评论
Guillaume
Guillaume 2020-2-4
But GD is 5x13 zero matrix.
Since nlat is 5 and nlon is 13, this is exactly what it should be. It's not clear why it is a problem. What size did you expect?

请先登录,再进行评论。


Guillaume
Guillaume 2020-2-5
编辑:Guillaume 2020-2-5
"GD should be the answer of pdf function. not all zero."
Ah, I thought you were complaining about the size of the output, not its content.
There's two problems about the code you've posted but first some advice:
Never end a if...elseif... without a plain else... to catch the conditions you may have forgotten, e.g.
if somecondition
%do something
elseif someothercondition
%do something else
elseif anothercondition
%do another thing
else
error('There''s a problem. Review the conditions');
end
If you'd done that you would have found that neither your if condition or your elseif condition was true. This explains why GD is all zeros, you never put anything in it.
Second advice, if your elseif condition is supposed to be the exact opposite of the if condition, then you don't need an elseif. You already know that if the if condition is not true, the opposite is always true. So instead of
if somecondition
%do something
elseif oppositeof_somecondition
%dosomething else
end
Just write
if somecondition
%do something
else
%do something else
end
This is clearer and crucially saves you from having to write the opposite condition correctly. That's your first problem, the opposite of A > 1 & C > 0 is not A < 1 & B < 0 but A <= 1 | B <= 0.
The second problem is that the expression
GD(i,j) = (1/((std(DIST{i,j})*sqrt(2*pi))))*exp(-((DIST{i,j}-mean(DIST{i,j})).^2)/(2*(std(DIST{i,j}).^2)))
can only work if DIST{i, j} is scalar. If it is not, then DIST{i,j}-mean(DIST{i,j}) will be a vector or matrix (same size as DIST{i, J}) and you'll be trying to store that into the scalar GD(i, j).
So you need to review your equation, it's still not right.
  3 个评论
Guillaume
Guillaume 2020-2-5
"I need to substract mean(DIST{i,j}) from every element of DIST{i,j}."
It looks like that's what your code is doing. If that is what you want, it is not possible to store that result as an element of a matrix since matrix elements must be scalar and the result of the above is going to be a vector or matrix the same size as DIST{i, j}.
Perhaps GD should be a cell array?

请先登录,再进行评论。

类别

Help CenterFile Exchange 中查找有关 Motion Modeling and Coordinate Systems 的更多信息

Community Treasure Hunt

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

Start Hunting!

Translated by