Effacer les filtres
Effacer les filtres

Can I smoothen the cdfplot() result?

8 vues (au cours des 30 derniers jours)
okoth ochola
okoth ochola le 11 Mar 2024
Commenté : Alexander le 26 Mar 2024
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 commentaires
Mathieu NOE
Mathieu NOE le 11 Mar 2024
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
Mathieu NOE le 25 Mar 2024
hello
so have you solved your problem ?

Connectez-vous pour commenter.

Réponse acceptée

Alexander
Alexander le 25 Mar 2024
@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 commentaires
okoth ochola
okoth ochola le 26 Mar 2024
Thank you all fro these wonderful answers. You dont how much have been helped
Alexander
Alexander le 26 Mar 2024
Cheers and upvote Mathieu also.

Connectez-vous pour commenter.

Plus de réponses (1)

Mathieu NOE
Mathieu NOE le 26 Mar 2024
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
end

Catégories

En savoir plus sur Get Started with Curve Fitting Toolbox dans Help Center et File Exchange

Produits


Version

R2018a

Community Treasure Hunt

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

Start Hunting!

Translated by