Curve fitting -twoway-

[pierre martin]

Hi,

I would like to construct a curve or mathematical function that has the best fit to a series of data points.

I am supposed to fit a curve of the form y=x^p through a plot of SD systolic blood pressure (y-axis) against mean systolic blood pressure (x-axis).

I first did a scatterplot of SD systolic blood pressure (y-axis) against mean systolic blood pressure (x-axis), for all individuals in my cohort. But I have some trouble to fit the curve as mentioned before.

How could I find this coefficient p?

I tried to use the - twoway qfit - command but I am not sure of the result.

Any advices?

Thank you very much for your valuable help,

Pierre

[Maarten Buis]

sysuse auto, clear
fp  <weight>, dim(1) : regress price <weight>
fp plot, res(res)

[Nick Cox]

I am more familiar with power laws or power functions of the form y = A x^b. These can be fitted in various ways and the best way to do it will depend on the structure of errors. But the most common method in my experience is to take logarithms

log y = log A + b log x

and use plain or vanilla regression on log y and log x. If we rewrite log y =: Y, log x =: X, log A =: a it can be seen that we are just fitting

Y = a + b X

Here's a silly example, absent an example of your own data.

. webuse grunfeld

. scatter invest mvalue

. gen ln_invest = ln(invest)

. gen ln_mvalue = ln(mvalue)

.. regress ln_invest ln_mvalue

      Source |       SS           df       MS      Number of obs   =       200
-------------+----------------------------------   F(1, 198)       =    532.11
       Model |  302.973353         1  302.973353   Prob > F        =    0.0000
    Residual |  112.736662       198  .569377083   R-squared       =    0.7288
-------------+----------------------------------   Adj R-squared   =    0.7274
       Total |  415.710015       199  2.08899505   Root MSE        =    .75457

------------------------------------------------------------------------------
   ln_invest |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
   ln_mvalue |   .9952899   .0431467    23.07   0.000     .9102039    1.080376
       _cons |   -2.13098   .2761548    -7.72   0.000    -2.675562   -1.586398
------------------------------------------------------------------------------

. scatter ln_invest ln_mvalue, ms(none) mla(company) mlabpos(0) || lfit ln_invest  ln_mvalue

Here the power b is 0.995: just fancy that. The prefactor A is

. di exp(-2.13098)
.11872089

and doesn't seem very interesting whenever differing units are involved.

Whether this is a good idea for blood pressure is hard to guess, but power laws make most sense when both variables vary over a few orders of magnitude, which I don't think is true for any kind of blood pressure.

twoway qfit fits a quadratic, which isn't in general the same thing.

[pierre martin]

Actually, I am studying visit to visit blood pressure variability (measures of variability for a set of n blood pressure measurements over time x1, x2 ..., xn) and I would like to calculate a variation independent of mean (VIM). To do that, I am supposed to fit a curve of the form y=kx^p through a plot of SD systolic blood pressure (y-axis) against mean systolic blood pressure (x-axis), for all individuals in my cohort.

The parameter p is estimated from the data and k is a constant which can be chosen such that the values of VIM are on the same scale as values of SD. For example, if M is the average value of mean systolic blood pressure in the cohort, then k = M^p and the value of VIM for any individual is given by VIM systolic blood pressure = (k SD / mean(x)^p).


Appendix_e-2 ROTHWELL.doc


I first did a scatterplot of SD systolic blood pressure (y-axis) against mean systolic blood pressure (x-axis), for all individuals in my cohort. But I have some trouble to fit the curve as mentioned before.

I don't know if this message helps to better explain what I am trying to do but I will read your responses.

Thank you for your help +++

Pierre

[Nick Cox]

Thanks for the detail. I confess I didn't focus on the detail on blood pressure. k SD / mean^p if that is hat you are seeking is a different beast. I would use nonlinear least squares for that.

[pierre martin]

Dear Nick,

I should have better explained it before. Yes this is what I am trying to dot. Exactly what is explained in the attached file Appendix ROTHWELL, their work has been published in the Lancet 2010 and I am trying to study blood pressure variability in a same way.

I have naver used before non linear least squares and I am a little bit confused with that since I have no real statistical background. I did a scatterplot of SD systolic blood pressure (y-axis) against mean systolic blood pressure (x-axis), for all individuals in my cohort. But I have some trouble to fit the curve as mentioned before.

Would you have some time to explain me a little bit more?

Thank you ever so much for your help,

Pierre

[Nick Cox]

I see no attachment here. You'll get better advice if you show us an example of the data. Please do read and study http://www.statalist.org/forums/help especially http://www.statalist.org/forums/help#stata

[pierre martin]

Yes of course. Here is an example of my dataset :

ctrpat is the patient number
visite if the visit number (0 to 7 visits)
idpat is the patient number with the visit number
A_SCMMSE is the cognitive performance score
A_AGE is the age
pasa is the systolic blood pressure.

* Example generated by -dataex-. To install: ssc install dataex
clear
input byte visite str9 ctrpat int pasa str10 idpat byte(A_SCMMSE A_AGE)
1 "DB0005001" 130 "DB00050011" 28 65
2 "DB0005001" 140 "DB00050012"  .  .
3 "DB0005001" 140 "DB00050013" 30  .
4 "DB0005001" 130 "DB00050014"  .  .
5 "DB0005001" 135 "DB00050015" 30  .
6 "DB0005001" 130 "DB00050016"  .  .
7 "DB0005001" 130 "DB00050017" 30  .
1 "DB0005002" 120 "DB00050021" 30 76
2 "DB0005002" 130 "DB00050022"  .  .
3 "DB0005002" 125 "DB00050023" 29  .
4 "DB0005002" 130 "DB00050024"  .  .
5 "DB0005002" 110 "DB00050025" 30  .
6 "DB0005002" 110 "DB00050026"  .  .
7 "DB0005002" 160 "DB00050027" 30  .
1 "DB0005003" 130 "DB00050031" 25 85
2 "DB0005003" 120 "DB00050032"  .  .
3 "DB0005003" 110 "DB00050033" 27  .
4 "DB0005003" 140 "DB00050034"  .  .
5 "DB0005003" 120 "DB00050035" 25  .
6 "DB0005003" 135 "DB00050036"  .  .
7 "DB0005003" 100 "DB00050037" 26  .
1 "DB0005004" 140 "DB00050041" 30 76
2 "DB0005004" 130 "DB00050042"  .  .
3 "DB0005004" 155 "DB00050043" 29  .
4 "DB0005004" 150 "DB00050044"  .  .
5 "DB0005004" 160 "DB00050045" 30  .
6 "DB0005004" 145 "DB00050046"  .  .
7 "DB0005004" 140 "DB00050047" 30  .
1 "DB0007001" 135 "DB00070011" 30 69
2 "DB0007001" 140 "DB00070012"  .  .
1 "DB0007002" 140 "DB00070021" 30 77
1 "DB0007003" 140 "DB00070031"  . 66
1 "DB0010001" 120 "DB00100011" 30 66
2 "DB0010001" 130 "DB00100012"  .  .
3 "DB0010001" 125 "DB00100013"  .  .
1 "DB0010002" 125 "DB00100021"  . 78
1 "DB0010003" 125 "DB00100031"  . 82
2 "DB0010003" 130 "DB00100032"  .  .
1 "DB0010004" 120 "DB00100041"  . 86
2 "DB0010004" 120 "DB00100042"  .  .
1 "DB0015001" 130 "DB00150011"  . 81
2 "DB0015001" 120 "DB00150012"  .  .
3 "DB0015001" 120 "DB00150013"  .  .
4 "DB0015001" 120 "DB00150014"  .  .
5 "DB0015001" 120 "DB00150015" 19  .
6 "DB0015001" 125 "DB00150016"  .  .
7 "DB0015001" 120 "DB00150017" 12  .
1 "DB0015002" 120 "DB00150021"  . 74
2 "DB0015002" 120 "DB00150022"  .  .
3 "DB0015002" 120 "DB00150023"  .  .
4 "DB0015002" 120 "DB00150024"  .  .
5 "DB0015002" 120 "DB00150025"  .  .
6 "DB0015002" 120 "DB00150026"  .  .
7 "DB0015002" 125 "DB00150027"  .  .
1 "DB0015003" 130 "DB00150031"  . 79
2 "DB0015003" 130 "DB00150032"  .  .
3 "DB0015003" 130 "DB00150033"  .  .
4 "DB0015003" 130 "DB00150034"  .  .
5 "DB0015003" 130 "DB00150035"  .  .
1 "DB0015004" 120 "DB00150041"  . 75
2 "DB0015004" 120 "DB00150042"  .  .
3 "DB0015004" 125 "DB00150043"  .  .
4 "DB0015004" 120 "DB00150044"  .  .
5 "DB0015004" 120 "DB00150045"  .  .
6 "DB0015004" 120 "DB00150046"  .  .
7 "DB0015004" 130 "DB00150047"  .  .
1 "DB0018001" 130 "DB00180011" 30 66
2 "DB0018001" 130 "DB00180012"  .  .
3 "DB0018001" 140 "DB00180013" 30  .
4 "DB0018001" 130 "DB00180014"  .  .
5 "DB0018001" 130 "DB00180015" 30  .
6 "DB0018001" 120 "DB00180016"  .  .
7 "DB0018001" 130 "DB00180017" 30  .
1 "DB0018002" 130 "DB00180021" 24 71
2 "DB0018002" 140 "DB00180022"  .  .
3 "DB0018002" 130 "DB00180023" 30  .
4 "DB0018002" 130 "DB00180024"  .  .
5 "DB0018002" 120 "DB00180025" 28  .
6 "DB0018002" 130 "DB00180026"  .  .
7 "DB0018002" 140 "DB00180027" 29  .
1 "DB0018003" 120 "DB00180031" 25 83
2 "DB0018003" 130 "DB00180032"  .  .
3 "DB0018003" 150 "DB00180033" 26  .
4 "DB0018003" 120 "DB00180034"  .  .
5 "DB0018003" 120 "DB00180035" 20  .
6 "DB0018003" 120 "DB00180036"  .  .
7 "DB0018003" 140 "DB00180037" 23  .
1 "DB0018004" 130 "DB00180041" 30 78
2 "DB0018004" 130 "DB00180042"  .  .
3 "DB0018004" 130 "DB00180043" 27  .
4 "DB0018004" 120 "DB00180044"  .  .
5 "DB0018004" 130 "DB00180045" 30  .
6 "DB0018004" 140 "DB00180046"  .  .
7 "DB0018004" 140 "DB00180047" 29  .
1 "DB0018005" 130 "DB00180051" 28 82
2 "DB0018005" 120 "DB00180052"  .  .
3 "DB0018005" 140 "DB00180053" 30  .
4 "DB0018005" 130 "DB00180054"  .  .
5 "DB0018005" 130 "DB00180055" 29  .
6 "DB0018005" 120 "DB00180056"  .  .
7 "DB0018005" 130 "DB00180057" 29  .
1 "DB0018006" 130 "DB00180061" 23 77
2 "DB0018006" 130 "DB00180062"  .  .
3 "DB0018006" 140 "DB00180063" 23  .
4 "DB0018006" 130 "DB00180064"  .  .
5 "DB0018006" 130 "DB00180065" 24  .
6 "DB0018006" 130 "DB00180066"  .  .
7 "DB0018006" 120 "DB00180067" 24  .
1 "DB0018007" 130 "DB00180071" 30 76
2 "DB0018007" 130 "DB00180072"  .  .
3 "DB0018007" 130 "DB00180073" 30  .
4 "DB0018007" 120 "DB00180074"  .  .
5 "DB0018007" 120 "DB00180075" 30  .
6 "DB0018007" 120 "DB00180076"  .  .
7 "DB0018007" 140 "DB00180077" 30  .
1 "DB0018008" 130 "DB00180081" 24 80
2 "DB0018008" 140 "DB00180082"  .  .
3 "DB0018008" 140 "DB00180083" 21  .
4 "DB0018008" 120 "DB00180084"  .  .
5 "DB0018008" 130 "DB00180085" 25  .
6 "DB0018008" 130 "DB00180086"  .  .
7 "DB0018008" 130 "DB00180087" 21  .
1 "DB0018009" 135 "DB00180091" 30 78
2 "DB0018009" 110 "DB00180092"  .  .
end

I am interested in studying the impact of visit to visit blood pressure variability over time on cognitive performances. I studied variability as a coefficient of variation, standard deviation and so on and I would like to study it as the "variation independent of mean".

Appendix_e-2 ROTHWELL.doc

Table e5 and paragraph "variation independent of mean"

This document (Rothwell Lancet 2010) explains how to calculate "variation independent of mean" as I tried to explain in the previous posts but I still have some trouble to make the calculation based on the scatterplot of SD systolic blood pressure (y-axis) against mean systolic blood pressure (x-axis).

Would you have some advices to help me?

Thank you ever so much for your help.

Pierre

[Nick Cox]

Sorry, but I won't undertake to read Word documents.

You cited VIM systolic blood pressure = (k SD / mean(x)^p).

which variables in your data are

VIM systolic blood pressure

SD

mean(x)

[pierre martin]

Of course I totally understand Nick. It was just in case my explanation was unclear.

VIM systolic blood pressure is the variable I want to construct.

SD is the variable "sdpasa" in my dataset.

mean(x) is the variable "meanpasa" in my dataset.


To do that, I am supposed to fit a curve of the form y=kx^p through a plot of SD systolic blood pressure (y-axis) against mean systolic blood pressure (x-axis), for all individuals in my cohort.

The parameter p is estimated from the data and k is a constant which can be chosen such that the values of VIM are on the same scale as values of SD. For example, if M is the average value of mean systolic blood pressure in the cohort, then k = M^p and the value of VIM for any individual is given by VIM systolic blood pressure = (k SD / mean(x)^p).

* Example generated by -dataex-. To install: ssc install dataex
clear
input byte visite long ctrpat int pasa str10 idpat byte(A_SCMMSE A_AGE) float(meanpasa sdpasa)
0 10005001 130 "DB00050011" 28 65 133.57143 4.7559485
1 10005001 140 "DB00050012"  .  . 133.57143 4.7559485
2 10005001 140 "DB00050013" 30  . 133.57143 4.7559485
3 10005001 130 "DB00050014"  .  . 133.57143 4.7559485
4 10005001 135 "DB00050015" 30  . 133.57143 4.7559485
5 10005001 130 "DB00050016"  .  . 133.57143 4.7559485
6 10005001 130 "DB00050017" 30  . 133.57143 4.7559485
0 10005002 120 "DB00050021" 30 76 126.42857 17.008402
1 10005002 130 "DB00050022"  .  . 126.42857 17.008402
2 10005002 125 "DB00050023" 29  . 126.42857 17.008402
3 10005002 130 "DB00050024"  .  . 126.42857 17.008402
4 10005002 110 "DB00050025" 30  . 126.42857 17.008402
5 10005002 110 "DB00050026"  .  . 126.42857 17.008402
6 10005002 160 "DB00050027" 30  . 126.42857 17.008402
0 10005003 130 "DB00050031" 25 85 122.14286 14.099983
1 10005003 120 "DB00050032"  .  . 122.14286 14.099983
2 10005003 110 "DB00050033" 27  . 122.14286 14.099983
3 10005003 140 "DB00050034"  .  . 122.14286 14.099983
4 10005003 120 "DB00050035" 25  . 122.14286 14.099983
5 10005003 135 "DB00050036"  .  . 122.14286 14.099983
6 10005003 100 "DB00050037" 26  . 122.14286 14.099983
0 10005004 140 "DB00050041" 30 76 145.71428 10.177005
1 10005004 130 "DB00050042"  .  . 145.71428 10.177005
2 10005004 155 "DB00050043" 29  . 145.71428 10.177005
3 10005004 150 "DB00050044"  .  . 145.71428 10.177005
4 10005004 160 "DB00050045" 30  . 145.71428 10.177005
5 10005004 145 "DB00050046"  .  . 145.71428 10.177005
6 10005004 140 "DB00050047" 30  . 145.71428 10.177005
0 10007001 135 "DB00070011" 30 69     137.5  3.535534
1 10007001 140 "DB00070012"  .  .     137.5  3.535534
0 10007002 140 "DB00070021" 30 77       140         .
0 10007003 140 "DB00070031"  . 66       140         .
0 10010001 120 "DB00100011" 30 66       125         5
1 10010001 130 "DB00100012"  .  .       125         5
2 10010001 125 "DB00100013"  .  .       125         5
0 10010002 125 "DB00100021"  . 78       125         .
0 10010003 125 "DB00100031"  . 82     127.5  3.535534
1 10010003 130 "DB00100032"  .  .     127.5  3.535534
0 10010004 120 "DB00100041"  . 86       120         0
1 10010004 120 "DB00100042"  .  .       120         0
0 10015001 130 "DB00150011"  . 81 122.14286  3.933979
1 10015001 120 "DB00150012"  .  . 122.14286  3.933979
2 10015001 120 "DB00150013"  .  . 122.14286  3.933979
3 10015001 120 "DB00150014"  .  . 122.14286  3.933979
4 10015001 120 "DB00150015" 19  . 122.14286  3.933979
5 10015001 125 "DB00150016"  .  . 122.14286  3.933979
6 10015001 120 "DB00150017" 12  . 122.14286  3.933979
0 10015002 120 "DB00150021"  . 74 120.71429 1.8898224
1 10015002 120 "DB00150022"  .  . 120.71429 1.8898224
2 10015002 120 "DB00150023"  .  . 120.71429 1.8898224
3 10015002 120 "DB00150024"  .  . 120.71429 1.8898224
4 10015002 120 "DB00150025"  .  . 120.71429 1.8898224
5 10015002 120 "DB00150026"  .  . 120.71429 1.8898224
6 10015002 125 "DB00150027"  .  . 120.71429 1.8898224
0 10015003 130 "DB00150031"  . 79       130         0
1 10015003 130 "DB00150032"  .  .       130         0
2 10015003 130 "DB00150033"  .  .       130         0
3 10015003 130 "DB00150034"  .  .       130         0
4 10015003 130 "DB00150035"  .  .       130         0
0 10015004 120 "DB00150041"  . 75 122.14286  3.933979
1 10015004 120 "DB00150042"  .  . 122.14286  3.933979
2 10015004 125 "DB00150043"  .  . 122.14286  3.933979
3 10015004 120 "DB00150044"  .  . 122.14286  3.933979
4 10015004 120 "DB00150045"  .  . 122.14286  3.933979
5 10015004 120 "DB00150046"  .  . 122.14286  3.933979
6 10015004 130 "DB00150047"  .  . 122.14286  3.933979
0 10018001 130 "DB00180011" 30 66       130  5.773503
1 10018001 130 "DB00180012"  .  .       130  5.773503
2 10018001 140 "DB00180013" 30  .       130  5.773503
3 10018001 130 "DB00180014"  .  .       130  5.773503
4 10018001 130 "DB00180015" 30  .       130  5.773503
5 10018001 120 "DB00180016"  .  .       130  5.773503
6 10018001 130 "DB00180017" 30  .       130  5.773503
0 10018002 130 "DB00180021" 24 71 131.42857  6.900656
1 10018002 140 "DB00180022"  .  . 131.42857  6.900656
2 10018002 130 "DB00180023" 30  . 131.42857  6.900656
3 10018002 130 "DB00180024"  .  . 131.42857  6.900656
4 10018002 120 "DB00180025" 28  . 131.42857  6.900656
5 10018002 130 "DB00180026"  .  . 131.42857  6.900656
6 10018002 140 "DB00180027" 29  . 131.42857  6.900656
0 10018003 120 "DB00180031" 25 83 128.57143 12.149858
1 10018003 130 "DB00180032"  .  . 128.57143 12.149858
2 10018003 150 "DB00180033" 26  . 128.57143 12.149858
3 10018003 120 "DB00180034"  .  . 128.57143 12.149858
4 10018003 120 "DB00180035" 20  . 128.57143 12.149858
5 10018003 120 "DB00180036"  .  . 128.57143 12.149858
6 10018003 140 "DB00180037" 23  . 128.57143 12.149858
0 10018004 130 "DB00180041" 30 78 131.42857  6.900656
1 10018004 130 "DB00180042"  .  . 131.42857  6.900656
2 10018004 130 "DB00180043" 27  . 131.42857  6.900656
3 10018004 120 "DB00180044"  .  . 131.42857  6.900656
4 10018004 130 "DB00180045" 30  . 131.42857  6.900656
5 10018004 140 "DB00180046"  .  . 131.42857  6.900656
6 10018004 140 "DB00180047" 29  . 131.42857  6.900656
0 10018005 130 "DB00180051" 28 82 128.57143  6.900656
1 10018005 120 "DB00180052"  .  . 128.57143  6.900656
2 10018005 140 "DB00180053" 30  . 128.57143  6.900656
3 10018005 130 "DB00180054"  .  . 128.57143  6.900656
4 10018005 130 "DB00180055" 29  . 128.57143  6.900656
5 10018005 120 "DB00180056"  .  . 128.57143  6.900656
6 10018005 130 "DB00180057" 29  . 128.57143  6.900656
0 10018006 130 "DB00180061" 23 77       130  5.773503
1 10018006 130 "DB00180062"  .  .       130  5.773503
2 10018006 140 "DB00180063" 23  .       130  5.773503
3 10018006 130 "DB00180064"  .  .       130  5.773503
4 10018006 130 "DB00180065" 24  .       130  5.773503
5 10018006 130 "DB00180066"  .  .       130  5.773503
6 10018006 120 "DB00180067" 24  .       130  5.773503
0 10018007 130 "DB00180071" 30 76 127.14286  7.559289
1 10018007 130 "DB00180072"  .  . 127.14286  7.559289
2 10018007 130 "DB00180073" 30  . 127.14286  7.559289
3 10018007 120 "DB00180074"  .  . 127.14286  7.559289
4 10018007 120 "DB00180075" 30  . 127.14286  7.559289
5 10018007 120 "DB00180076"  .  . 127.14286  7.559289
6 10018007 140 "DB00180077" 30  . 127.14286  7.559289
0 10018008 130 "DB00180081" 24 80 131.42857  6.900656
1 10018008 140 "DB00180082"  .  . 131.42857  6.900656
2 10018008 140 "DB00180083" 21  . 131.42857  6.900656
3 10018008 120 "DB00180084"  .  . 131.42857  6.900656
4 10018008 130 "DB00180085" 25  . 131.42857  6.900656
5 10018008 130 "DB00180086"  .  . 131.42857  6.900656
6 10018008 130 "DB00180087" 21  . 131.42857  6.900656
0 10018009 135 "DB00180091" 30 78 125.71429  8.380817
1 10018009 110 "DB00180092"  .  . 125.71429  8.380817
end

Thank you so much for your help,

Pierre

[Nick Cox]

OK, If I understand correctly you want something like

egen tag = tag(ctrpat) 

gen ln_sd = ln(sd) 
gen ln_mean = ln(mean) 

regress ln_sd ln_mean if tag 
predict predict 
scatter ln_sd ln_mean if tag, ytitle(ln SD) yla(, ang(h)) ms(Oh) || ///
line predict ln_mean if tag, sort xtitle(ln mean) aspect(1) legend(off)

[pierre martin]

Nick,

It is perfectly what I was looking for and I created my variation independent of mean variable. It helped me to better deal with logarithmic transformations.

Thank you so much again for the time you allowed me and your help.

Have a nice day,

Best,

Pierre