Can I smoothen the cdfplot() result?

Question

okoth ochola 2024-3-11

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2092911-can-i-smoothen-the-cdfplot-result

评论： Alexander 2024-3-26

Hi, suppose I have a sample data x, and i want to plot the the cdf function on the data and plot it, ofcoz I will use cdfplot(x) and I will a graph similar to the one attached below. My question is, can i smoothen it, is there a way I can make the plot smoother and not look like a a step funtion plot?

3 个评论
显示 1更早的评论隐藏 1更早的评论

Mathieu NOE 2024-3-11

even with smoothing , I doubt you can make those large flats disappear

wuld probably rather try to fit either a polynomial or exponential function

as mentionned by @Alexander, sharing the data would be a plus...

Mathieu NOE 2024-3-25

hello

so have you solved your problem ?

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Alexander 2024-3-25

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2092911-can-i-smoothen-the-cdfplot-result#answer_1430826

@Mathieu NOE, thank you for remaindig me. This question kept nagging me, but I just forgotten it. Here are some rudimentary aproches. As the data were artificially generated, with real measurement the solution are only a rough direction what can be done and needed a lot of refinement.

commandwindow;

y = [];

for ii = 1:20

y = [y ones(1,100*ii)*ii];

end

n = length(y);

x = linspace(1,25,n);

plot(x,y);grid minor;

% 1. Polynomial Fit, 2. Order

% As suggested by @Mathieu NOE

P = polyfit(x,y,2);

yNew = polyval(P,x);

plot(x,y,x,yNew);grid minor; title('Polyfit');pause;

% 2. Logarithmic Fit

% Could be done much easier with lsqcurvefit, if you have optimization

% toobox. I don't have it, hence, it's handcrafted.

yLogSum = sum(y.*log(x)); ySum = sum(y); sumLnx = sum(log(x));

Zaeler = yLogSum - (ySum*sumLnx)/n;

Nenner = sum(log(x).^2) - sumLnx^2/n;

b = Zaeler/Nenner;

a = ySum/n - b*sumLnx/n;

yLogRes = a + b*log(x);

plot(x,y,x,yLogRes);grid minor; title('Logarithmic Fit');pause;

% 3. Filtering

% I just used a butterworth filter. As more longer the steps are, as higher

% is the weaviness of the output. The filter parameters has to be adjusted

% to fit the input curve. Another filter type might fit better.

[B,A] = butter(2,0.0005);

yFF = filtfilt(B,A,y);

plot(x,y,x,yFF);grid minor; title('Butterworth filtered');pause;

% 4. Step/Edge Evaluation

% The input was generated. With measured data the following lines must be

% adjusted.

% Also the data has to be interpolated between the x-values. yDiff has now

% as many data values as there are steps in the measurement. Maybe with

% interp().

yDiff = diff(y);

yFind = find(yDiff > 0.5);

Edge = 1; %Upper Edge = 1; lower edge = 0;

plot(x,y,x(yFind+Edge),y(yFind+Edge));grid minor; title('Edge Evaluation');

legend('Original',['Edge ' num2str(Edge)], 'location','northwest')

4 个评论
显示 2更早的评论隐藏 2更早的评论

Mathieu NOE 2024-3-26

在 MATLAB Online 中打开

hello again

I slightly modified your code (first I reduced the data size by factor 10 to run it faster, but that's really minor) , then I think using the upper and lower edge approach we can build a good solution

would like to continue with real data

y = [];

for ii = 1:20

y = [y ones(1,10*ii)*ii];

end

n = length(y);

x = linspace(1,25,n);

% 1. Polynomial Fit, 2. Order

% As suggested by @Mathieu NOE

P = polyfit(x,y,6);

yNew = polyval(P,x);

figure,plot(x,y,x,yNew);grid minor; title('Polyfit');

% 2. Logarithmic Fit

% Could be done much easier with lsqcurvefit, if you have optimization

% toobox. I don't have it, hence, it's handcrafted.

yLogSum = sum(y.*log(x)); ySum = sum(y); sumLnx = sum(log(x));

Zaeler = yLogSum - (ySum*sumLnx)/n;

Nenner = sum(log(x).^2) - sumLnx^2/n;

b = Zaeler/Nenner;

a = ySum/n - b*sumLnx/n;

yLogRes = a + b*log(x);

figure,plot(x,y,x,yLogRes);grid minor; title('Logarithmic Fit');

% 3. Filtering

% I just used a butterworth filter. As more longer the steps are, as higher

% is the weaviness of the output. The filter parameters has to be adjusted

% to fit the input curve. Another filter type might fit better.

[B,A] = butter(2,0.005);

yFF = filtfilt(B,A,y);

plot(x,y,x,yFF);grid minor; title('Butterworth filtered');

% 4. Step/Edge Evaluation

% The input was generated. With measured data the following lines must be

% adjusted.

% Also the data has to be interpolated between the x-values. yDiff has now

% as many data values as there are steps in the measurement. Maybe with

% interp().

yDiff = diff(y);

yFind = find(yDiff > 0.5);

Edge = 0; %Upper Edge = 1; lower edge = 0;

xl = x(yFind+Edge);

yl = y(yFind+Edge);

Edge = 1; %Upper Edge = 1; lower edge = 0;

xu = x(yFind+Edge);

yu = y(yFind+Edge);

% Median curve

xm = (xu+xl)/2;

ym = (yu+yl)/2;

figure,plot(x,y,xl,yl,xu,yu,xm,ym);grid minor; title('Edge Evaluation');

legend('Original','Lower Edge','Upper Edge','Median curve' , 'location','northwest')

okoth ochola 2024-3-26

Thank you all fro these wonderful answers. You dont how much have been helped

Alexander 2024-3-26

Cheers and upvote Mathieu also.

请先登录，再进行评论。

Answer 2

Mathieu NOE 2024-3-26

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2092911-can-i-smoothen-the-cdfplot-result#answer_1431161

在 MATLAB Online 中打开

some more code to look at - I kept only the best solutions , even though there are probably lot of other solutions

also I introduce some randomness in the data to see how robust is the code

enjoy !

y = [];

for ii = 1:20

y = [y ones(1,round(10*ii+25*rand))*ii];

end

n = length(y);

x = linspace(1,25,n);

% add some noise on y

y = y + 0.05*randn(size(y));

% 1. Polynomial Fit,

% As suggested by @Mathieu NOE

P = polyfit(x,y,6);

yNew = polyval(P,x);

figure,plot(x,y,x,yNew);grid minor; title('Polyfit');

% 2. Step/Edge Evaluation

% The input was generated. With measured data the following lines must be

% adjusted.

% Also the data has to be interpolated between the x-values. yDiff has now

% as many data values as there are steps in the measurement. Maybe with

% interp().

yDiff = diff(y);

yFind = find(yDiff > max(yDiff)/2);

Edge = 0; %Upper Edge = 1; lower edge = 0;

xl = x(yFind+Edge);

yl = y(yFind+Edge);

Edge = 1; %Upper Edge = 1; lower edge = 0;

xu = x(yFind+Edge);

yu = y(yFind+Edge);

% Median curve

xm = (xu+xl)/2;

ym = (yu+yl)/2;

figure,plot(x,y,xl,yl,xu,yu,xm,ym);grid minor; title('Edge Evaluation');

legend('Original','Lower Edge','Upper Edge','Median curve' , 'location','northwest')

% 3. Envelop (secant method)

N = 200;

[yu] = env_secant(x, y, N, 'top'); % (see function attached)

[yl] = env_secant(x, y, N, 'bottom'); % (see function attached)

% Median curve

ym = (yu+yl)/2;

figure,plot(x,y,x,yl,x,yu,x,ym);grid minor; title('Envelop (secant method)');

legend('Original','Lower Envelop','Upper Envelop','Median curve' , 'location','northwest')

%%%%%%%%%%%%%%%%%%%%%%%%%%%

function [env] = env_secant(x_data, y_data, view, side)

% Function call: env_secant(x_data, y_data, view, side)

% Calculates the top envelope of data <y_data> over <x_data>.

% Method used: 'secant-method'

% env_secant() observates the max. slope of about <view> points,

% and joints them to the resulting envelope.

% An interpolation over original x-values is done finally.

% <side> ('top' or 'bottom') defines which side to evolve.

% Author: Andreas Martin, Volkswagen AG, Germany

if nargin == 0

test( false );

test( true );

return

end

if nargin < 4

error( '%s needs at least 4 input arguments!', mfilename );

end

assert( isnumeric(view) && isscalar(view) && view > 1, ...

'Parameter <view> must be a value greater than 1!' );

assert( isvector(x_data) && isnumeric(x_data) && all( isfinite( x_data ) ), ...

'Parameter <x_data> has to be of vector type, holding finite numeric values!' );

assert( isvector(y_data) && isnumeric(y_data) && all( isfinite( y_data ) ), ...

'Parameter <y_data> has to be of vector type, holding finite numeric values!' );

assert( isequal( size(x_data), size(y_data) ), ...

'Parameters <x_data> and <y_data> must have same size and dimension!' );

assert( ischar(side), ...

'Parameter <side> must be ''top'' or ''bottom''!' );

switch lower( side )

case 'top'

side = 1;

case 'bottom'

side = -1;

otherwise

error( 'Parameter <side> must be ''top'' or ''bottom''!' );

end

sz = size( x_data );

x_data = x_data(:);

x_diff = diff( x_data );

x_diff = [min(x_diff), max(x_diff)];

assert( x_diff(1) > 0, '<x_data> must be monotonic increasing!' );

y_data = y_data(:);

data_len = length( y_data );

x_new = [];

y_new = [];

if diff( x_diff ) < eps( max(x_data) ) + eps( min(x_data) )

% x_data is equidistant

search_fcn = @( y_data, ii, i ) ...

max( ( y_data(ii) - y_data(i) ) ./ (ii-i)' * side );

else

% x_data is not equidistant

search_fcn = @( y_data, ii, i ) ...

max( ( y_data(ii) - y_data(i) ) ./ ( x_data(ii) - x_data(i) ) * side );

end

i = 1;

while i < data_len;

ii = i+1:min( i + view, data_len );

[ m, idx ] = search_fcn( y_data, ii, i );

% New max. slope: store new "observation point"

i = i + idx;

x_new(end+1) = x_data(i);

y_new(end+1) = y_data(i);

end;

env = interp1( x_new, y_new, x_data, 'linear', 'extrap' );

env = reshape( env, sz );

function test( flagMonotonic )

npts = 100000;

y_data = cumsum( randn( npts, 1 ) ) .* cos( (1:npts)/50 )' + 100 * cos( (1:npts)/6000 )';

if flagMonotonic

x_data = (1:npts)' + 10;

else

x_diff = rand( size( y_data ) );

x_data = cumsum( x_diff );

end

view = ceil( npts * 0.01 ); % 1 Percent of total length

env_up = env_secant( x_data, y_data, view, 'top' );

env_lo = env_secant( x_data, y_data, view, 'bottom' );

figure

plot( x_data, y_data, '-', 'Color', 0.8 * ones(3,1) );

hold all

h(1) = plot( x_data, env_up, 'b-', 'DisplayName', 'top' );

h(2) = plot( x_data, env_lo, 'g-', 'DisplayName', 'bottom' );

grid

legend( 'show', h )

end

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

Can I smoothen the cdfplot() result?

3 个评论
显示 1更早的评论隐藏 1更早的评论

采纳的回答

4 个评论
显示 2更早的评论隐藏 2更早的评论

更多回答（1 个）

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

Can I smoothen the cdfplot() result?

3 个评论 显示 1更早的评论隐藏 1更早的评论

采纳的回答

4 个评论 显示 2更早的评论隐藏 2更早的评论

更多回答（1 个）

0 个评论 显示 -2更早的评论隐藏 -2更早的评论

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

3 个评论
显示 1更早的评论隐藏 1更早的评论

4 个评论
显示 2更早的评论隐藏 2更早的评论

0 个评论
显示 -2更早的评论隐藏 -2更早的评论