Finding Steepest Linear Portion of Data

I've been trying to use fitlm as a way of finding the steepest linear portion of a graph but cannot get it to work correctly. Perhaps there is a better way I could be finding this out?

Essentially I have data (stress and strain) which looks like the following:

I'm trying to use Matlab to automatically find the steepest linear portion of the graph. This is what I have been trying so far using the fitlm function.

    R=0;

for x = 100:size(strain)

      fitlm(strain(10:x),stress(10:x));

      if ans.Rsquared.Ordinary > R
          R=ans.Rsquared.Ordinary;
          index=x;
      end
end

I can get the code to go along the data and find out where the r-squared value seems to fit best and identify a linear region of the data, but how can I get Matlab to identify the steepest linear region?

This is the output of plotting from the start of the data to the index point:

[Oliver Ramos]

Hi Jonathan,

I have the same problem, I need to find the linear portion from my data points (quadratic) by finding the best R^2 based on a fixed window (number of data that forms a line). May I know how did you plot the resulting best linear portion from the code?

R=0;

for x = 100:size(strain)

fitlm(strain(10:x),stress(10:x));

if ans.Rsquared.Ordinary > R

R=ans.Rsquared.Ordinary;

index=x;

end

end

Answer 1

[Star Strider]

You have a non-linear model, so fitting linear segments is likely not going to be easy. The best way to fit your data is to use a model of the process that created your data (if you have an analytical version of it) and then a non-linear regression, such as fitnlm or nlinfit. Then take the analytic derivative of the fitted function.

If you’re looking for linear segments in your data and don’t want to fit it, I would start by using the gradient function to calculate the numerical derivative of your data (since they appear reasonably noise-free). This should produce a series of ‘steps’ that are flat in the linear sections. Choose the one that meets your criteria in terms of length and ‘flatness’, and use it. One of the histogram functions could be helpful here, since you can define the ranges of the bins. This will eliminate a lot of tedious threshold programming. (I don’t have your data, so I can’t write code to analyse it.)

[Jonathan Acheson]

Thank you for your help. A work colleague also suggested using a window size to calculate gradients and then find the highest. The following code is what I used to achieve this:

load('strain.mat')
load('stress.mat')
%  n =  number of points between starting point and ending point
n=150;

grad=zeros(585-(n-1),1);

for i=1:585-(n-1)

      grad(i,1)=(stress(i+(n-1))-stress(i))/(strain(i+(n-1))-strain(i));

end

ind=find(grad==max(grad));

x_strain=[strain(ind),strain(ind+(n-1))];

y_stress=[stress(ind),stress(ind+(n-1))];

figure

plot(strain,stress,'r')

hold on

plot(x_strain,y_stress,'b')

c=y_stress(1)-(grad(ind)*x_strain(1));

x_intercept=-c/grad(ind);

hold on

x_int=[x_intercept,strain(ind+(n-1))];

y_int=[0,stress(ind+(n-1))];

plot(x_int,y_int,'g')

With the resulting graph looking like:

With the green line being the steepest linear region, I just need to keep an eye on the window size.

[Image Analyst]

Well that's what I suggested but I think my method was more robust since your colleague's method just looks at two adjacent elements to compute the gradient, "grad". This can give very high slopes in some areas, like around a strain of .03 and .04 where it is very very steep but just for an very very short distance. With my method you can have control over the window width - it's not just 1 element.

[Star Strider]

My pleasure.

[Abdulaziz Abutunis]

Hi Star Strider, I have the same problem, I have a set of non-linear data that has a linear segment close to the lift end of the curve. I wonder how would you suggest to use the gradient function to locate the best range of data that represents the linear segment and thin find the slop and the interception point. The data are as follow and the linear segment always crosses y coordinate means that the linear segment should have negative and positive elements from x array and usually between -4 and 5 but sometimes this range slightly increases or decreases.

X=[ -5 -4 -3 -2 -1 0 1 2 3 4 5 6 7 8 9 10];

y=[ -0.2757 -0.3914 -0.4219 -0.3808 -0.3122 -0.1272 0.0800 0.1745 0.2545 0.3270 0.3932 0.4539 0.5091 0.5590 0.6040 0.6450];

Thanks in advance Aziz

Answer 2

[Image Analyst]

Steepest over what range of x? From one element to the next? Over 5 or 10 elements? Can you attach numerical data?

How about just scanning the data and fitting a line over the range, and keeping track of which location has the highest slope? Something like (untested):

maxSlope = -inf;
windowWidth = 5; % Whatever...
maxk = 1;
for k = 1 : length(x) - windowWidth 
    % Get a subset of the data that is just within the window.
    windowedx = x(k:k+windowWidth-1);
    windowedy = y(k:k+windowWidth-1);
    % Fit a line to that chunk of data.
    coefficients = polyfit(windowedx, windowedy, 1);
    % Get the slope.
    slope = coefficients(1);
    % See if this slope is steeper than any prior slope.
    if slope > maxSlope
        % It is steeper.  Keep track of this one.
        maxSlope = slope
        maxk = k;
    end
end