topic Re: Calculating area between two spline fits (bivariate) in Discussions

Calculating area between two spline fits (bivariate)

ms — Tue, 19 Jan 2010 00:03:14 GMT

I have used fit spline of the bivariate platform to parameterize fluorometric spectral curves and saved the prediction formula for my samples as well as blanks. I want to estimate the area under the peak, or more precisely, the area between the curves of the sample and the blank within a specified wavelength interval.

I could certainly do this in matlab with some effort, but it would be neat keep to the entire workflow within JMP. I think my approach would involve summing up slices which geometry is determined by the spline prediction formula matrix. However I am not sure how to interpret and use the prediction formula.

Any ideas that may lead me in the right direction on how to tackle this problem using JMP/JSL? Thanks!

Re: Calculating area between two spline fits (bivariate)

mpb — Wed, 20 Jan 2010 22:15:40 GMT

From the 8.0.1 Stat and Graph Guide:

"The cubic spline method uses a set of third-degree polynomials spliced together such that the resulting curve is continuous and smooth at the splices (knot points). The estimation is done by minimizing an objective function that is a combination of the sum of squares error and a penalty for curvature integrated over the curve extent. See the paper by Reinsch (1967) or the text by Eubank (1988) for a description of this method."

When you fit a spline in Fit Y by X there is a red triangle pull down menu with an item to Save Coefficients. This action creates a new table with columns X, A, B, C, D. Column X contains the unique X values from the data in ascending sort order. The other 4 columns are the coefficients of the third-degree polynomial used in an interval starting with the value of X for the same row up to the value of X for the next row. If the value for X in the same row is X1 then the polynomial for the interval is:

A + B*(X-X1) + C*(X-X1)^2 + D*(X-X1)^3

which is used for X up to X2, the next "knot". To find the exact area under this curve between x1 and x2, simply integrate this polynomial from X1 to X2 to get:

A*(X2-X1) + (B*(X2-X1)^2)/2 + (C*(X2-X1)^3)/3 + (D*(X2-X1)^4)/4

You need to make a new column called Area Slice and put in this formula. In the JMP formula, use X for X1 and Lag(X,-1) for X2. Then you will have a column of Area Slices. You can then add all of the area slices to get the total area under the curve from start to finish. You can also make a column of cumulative areas if you want to.

This gives you an exact value, not an approximation, for the area under the spline.

This can probably be scripted but I haven't tried.

To get the area between 2 such curves, simply find the area this way for both and subtract one from the other.

Michael

Re: Calculating area between two spline fits (bivariate)

ms — Fri, 22 Jan 2010 08:52:26 GMT

Thanks! It to works perfectly and am now trying to do it in JSL. I can post the script here when done, if anyone's interested.

Btw, there is a problem with the cumulative column. If do it manually by adding the formula
If(Row() == 1, :Slice Area, :Cumulative Area[Row() - 1] + :Slice Area[Row()]) to a column "Cumulative Area" it works. However if I do it by scripting using the command << set formula() I get an error "Illegal Reference (Recursion) at row....". The column is created but without the formula.

Re: Calculating area between two spline fits (bivariate)

mpb — Fri, 22 Jan 2010 12:45:39 GMT

Search this forum using: Column Formula Variable or something similar to find an example of how to successfully create a column formula via scripting which references an unknown-in-advance column. I think that will help.

Re: Calculating area between two spline fits (bivariate)

mpb — Sun, 24 Jan 2010 18:49:52 GMT

On JMP 8.0.2 under WinXP, using a sample JMP generated Spline coefficient table when prompted, the following script did work for me without error:

Re: Calculating area between two spline fits (bivariate)

ms — Mon, 25 Jan 2010 22:46:42 GMT

This script works, which made wonder why mine did not, as it is supposed to be identical to yours. When l compared them I found a typo in mine. Sorry about the confusion. And thanks for your help! I have a now working script for batch processing my fluorometry data.

The only annoying thing remaining is the Maximize() command that breaks the script if it fails to converge, i.e. find a local maxima. I have not found out how to interact with the error output in JSL (Like "if result==error, then change tolerance and repeat, else Max=result" or similar).

Thanks again!