1
$\begingroup$

I have a set of data that I think forms a power law relationship, however I am struggling to work out the equation of the relationship.

Here is a subset of the data I am working with:

df <- data.frame('x' = c(0.25, 1, 2, 2.75, 4.25, 6.5, 8, 13.25, 16.25, 19.25, 19.75, 26.5, 31, 37.75),
                 'y' = c(4.485605e-08, 1.430240e-08, 7.638950e-09, 6.776308e-09, 3.269885e-09, 
                     2.609455e-09, 4.378785e-09, 2.260540e-09, 2.039074e-09, 7.119317e-10,
                     2.252598e-09, 1.617082e-09, 7.511261e-09, 1.519275e-09))

Ploting this, I thought that the graph resembles that of y = 1/x:

enter image description here

I then plotted the log-log of the graph to see if it produced a straight line (indicating a power relationship) and attempted to fit a linear model to it to find the slope and the intercept.

enter image description here

The slope of the line was found to be -0.6099 and the intercept -18.178. However, when I put this back into the equation to work out the value of y for a given value of x it is way off.

For example if x=0.25 then y = (10^-18.178)*(0.25^-0.601).

This calculates y to be 1.527x10^-18 when the value of y for that data point was actually 4.486x10^-8.

I thought that perhaps is it simply that the heteroscedacity of the varience violates an assumption of the linear regression such that I can't use it for anything more than calculating the slope.

Am I correct, or is there something I've missed? Why is this value so far off?

$\endgroup$
  • 1
  • 1
    $\begingroup$ Thank you, I have seen that one links a paper explaining better ways to calculate the coefficients for power law relationships which I will check out. $\endgroup$ – tom91 Feb 14 '19 at 11:50
  • 2
    $\begingroup$ You're fitting using natural logarithms but using a power of 10 in the intercept calculations. $\endgroup$ – Nick Cox Feb 14 '19 at 11:50
  • $\begingroup$ @NickCox I have literally just realised this after working through yoav_aaa's answer below and now feel rather foolish. $\endgroup$ – tom91 Feb 14 '19 at 11:53
  • $\begingroup$ Happens to us all. Ideally you'd now just delete the thread. But there is an upvoted accepted answer (which doesn't really explain the issue!). $\endgroup$ – Nick Cox Feb 14 '19 at 11:58
1
$\begingroup$

I did the following calculation and got a different y value for x=0.25.
Intercept and slope are similar to OP's question.

x = (0.25, 1, 2, 2.75, 4.25, 6.5, 8, 13.25, 16.25, 19.25, 19.75, 26.5, 31, 37.75)
y = (4.485605e-08, 1.430240e-08, 7.638950e-09, 6.776308e-09, 3.269885e-09, 
                 2.609455e-09, 4.378785e-09, 2.260540e-09, 2.039074e-09, 7.119317e-10,
                 2.252598e-09, 1.617082e-09, 7.511261e-09, 1.519275e-09)
my_df = pd.DataFrame({'x':x , 'y': y})
my_df.loc[:, 'log_x'] = map(lambda x: math.log(x), my_df['x'])
my_df.loc[:, 'log_y'] = map(lambda x: math.log(x), my_df['y'])
model = LinearRegression()
model.fit(my_df['log_x'].values.reshape(-1, 1), my_df['log_y'].values)
model.coef_
[-0.60989947]
model.intercept_
[-18.17749636096132]
X = 0.25
log_y = model.intercept_ + model.coef_ * math.log(X)
y = math.exp(log_y)
print(y)
2.970364e-08
$\endgroup$
  • $\begingroup$ Interesting, could you post the values of model.coef_ and model.intercept_ please? I want to see how they vary from those calculate in the model I produced in r $\endgroup$ – tom91 Feb 14 '19 at 11:45
  • $\begingroup$ I have just realised the cause of the issue was me using a power of 10 after using natural logs. D'OH! $\endgroup$ – tom91 Feb 14 '19 at 11:56

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.