commons-user mailing list archives

Site index · List index
Message view « Date » · « Thread »
Top « Date » · « Thread »
From SUJIT PAL <sujit....@comcast.net>
Subject Re: [math] Getting second derivative for discrete dataset
Date Sat, 31 Mar 2012 17:59:14 GMT
Thank you Luc, this is very useful. Definitely much easier than what I was thinking of doing
:-).

-sujit

On Mar 31, 2012, at 10:13 AM, Luc Maisonobe wrote:

> Le 31/03/2012 01:31, SUJIT PAL a écrit :
>> Thanks Luc, thanks for getting back so quickly, and will wait for your information
tomorrow.
> 
> Hi Sujit,
> 
> Sorry for the delay.
> 
> If you know your sample points (xi, yi) are all separated by the same
> interval h (i.e x(i+1) = x(i) + h for all i), then you can compute any
> derivatives using the following method:
> 
> First write down the series expansion of the function y = f(x):
> 
>  y(x + k * h) = y(x)
>               + k * h * y'(x)
>               + (k * h)^2 / 2 y''(x)
>               + ...
>               + (k * h)^n / n! yn(x)
>               + ...
> 
> Then you truncate this expression to at least include the derivative you
> want, and apply it to a number of points equal to your truncation order.
> Here, as an example, we will truncate just to order 2 in order to have
> y'', so we have three components (y, y' and y''), so we need three
> points. Let's take k = -1, k = 0 and k = 2. We have the following
> expressions:
> 
> y(x - h) = y(x) - h y'(x) + h^2/2 y''(x)
> y(x)     = y(x)
> y(x+h)   = y(x) + h y'(x) + h^2/2 y''(x)
> 
> Note that the first and last expression differ only in the sign of the
> odd terms h y'(x), the even terms (h^0 and h^2) have the same sign.
> 
> You can rewrite these expressions as a simple square matrix vector
> multiplication:
> 
> [ y(x-h) ]   [ 1  -h  h^2/2 ]   [  y(x)  ]
> [ y(x)   ] = [ 1   0  0     ] * [  y'(x) ]
> [ y(x+h) ]   [ 1   h  h^2/2 ]   [ y''(x) ]
> 
> Then you invert this system and get:
> 
> [  y(x)  ]    [    0,        1,        0    ]   [ y(x-h) ]
> [  y'(x) ]  = [ -1 / 2h,     0,      1 / 2h ] * [ y(x)   ]
> [ y''(x) ]    [  1 / h^2, -2 / h^2,  1 / h^2]   [ y(x+h) ]
> 
> The third row is the one you want:
> 
> y''(x) = (y(x-h) - 2 * y(x) + y(x+h)) / h^2
> 
> You can check this formula using again the initial expansion, and you
> will see this combination basically canels out the y(x) and y'(x) part
> and only y''(x) remains.
> 
> You could also have used more points, which could improve accuracy as it
> fits your function to an higher degree. Beware that you should *not*
> increase the number of points too much, otherwise you would get Gibbs
> oscillation problems (i.e a polynomial that has very large changes
> between the sampling points).
> 
> Hope this helps
> Luc
> 
> 
>> 
>> -sujit
>> 
>> On Mar 30, 2012, at 2:36 PM, Luc Maisonobe wrote:
>> 
>>> 
>>> 
>>> 
>>> SUJIT PAL <sujit.pal@comcast.net> a écrit :
>>> 
>>>> Hi,
>>>> 
>>>> I have a (newbie) question about how to go about solving the problem
>>>> below with commons-math.
>>>> 
>>>> 1) I have histogram data (equal x intervals) based off some
>>>> distribution.
>>>> 2) I need a way to calculate the second differential between two given
>>>> x points in the distribution.
>>>> 
>>>> Here is what I think I should do:
>>>> 1) Feed the (x,y) values represented by the top of each histogram entry
>>>> within my x-range into a UnivariateRealInterpolator between ranges
>>>> (x1,y1) and (x2,y2).
>>>> 2) The interpolate method gives me back a PolynomialSplineFunction
>>>> (call if psf).
>>>> 3) I get the second differential UnivariateRealFunction (call it urf)
>>>> by calling psf.getPolynomialSplineDerivative().derivative()
>>>> 4) I compute the second derivative value as the difference of
>>>> urf.value(x2) - urf.value(x1)
>>>> 
>>>> Does this sound like a reasonable/accurate approach? Any
>>>> suggestions/gotchas or pointers to better approaches?
>>> 
>>> No, this is not an approach I would recommend. A much more straightforward
>>> way to compute this derivative is to apply a dedicated finite differences formula
>>> on three consecutive points. I'll give ou such a formula tomorrow (I am on a
mobile
>>> device right now and it is difficult to write a complète mail...)
>>> 
>>> Luc
>>> 
>>>> 
>>>> Thanks very much,
>>>> 
>>>> Sujit
>>>> 
>>>> 
>>>> ---------------------------------------------------------------------
>>>> To unsubscribe, e-mail: user-unsubscribe@commons.apache.org
>>>> For additional commands, e-mail: user-help@commons.apache.org
>>> 
>>> 
>>> ---------------------------------------------------------------------
>>> To unsubscribe, e-mail: user-unsubscribe@commons.apache.org
>>> For additional commands, e-mail: user-help@commons.apache.org
>>> 
>> 
>> 
>> ---------------------------------------------------------------------
>> To unsubscribe, e-mail: user-unsubscribe@commons.apache.org
>> For additional commands, e-mail: user-help@commons.apache.org
>> 
> 
> 
> ---------------------------------------------------------------------
> To unsubscribe, e-mail: user-unsubscribe@commons.apache.org
> For additional commands, e-mail: user-help@commons.apache.org
> 


---------------------------------------------------------------------
To unsubscribe, e-mail: user-unsubscribe@commons.apache.org
For additional commands, e-mail: user-help@commons.apache.org


Mime
View raw message