fitSuper: Data-Fitting Function for the Superellipse Equation

Description

fitSuper is used to estimate the parameters of the superellipse equation.

Usage

fitSuper(x, y, ini.val, control = list(), par.list = FALSE, 
         stand.fig = TRUE, angle = NULL, fig.opt = FALSE, np = 2000,
         xlim = NULL, ylim = NULL, unit = NULL, main = NULL)

Value

par: the estimated values of the parameters in the superellipse equation.
scan.length: the observed length of the polygon.
scan.width: the observed width of the polygon.
scan.area: the observed area of the polygon.
r.sq: the coefficient of determination between the observed and predicted polar radii.
RSS: the residual sum of squares between the observed and predicted polar radii.
sample.size: the number of data points used in the data fitting.
phi.stand.obs: the polar angles at the standard state.
phi.trans: the transferred polar angles rotated as defined by the user.
r.stand.obs: the observed polar radii at the standard state.
r.stand.pred: the predicted polar radii at the standard state.
x.stand.obs: the observed $x$ coordinates at the standard state.
x.stand.pred: the predicted $x$ coordinates at the standard state.
y.stand.obs: the observed $y$ coordinates at the standard state.
y.stand.pred: the predicted $y$ coordinates at the standard state.
r.obs: the observed polar radii at the transferred polar angles as defined by the user.
r.pred: the predicted polar radii at the transferred polar angles as defined by the user.
x.obs: the observed $x$ coordinates at the transferred polar angles as defined by the user.
x.pred: the predicted $x$ coordinates at the transferred polar angles as defined by the user.
y.obs: the observed $y$ coordinates at the transferred polar angles as defined by the user.
y.pred: the predicted $y$ coordinates at the transferred polar angles as defined by the user.

Arguments

x: the $x$ coordinates of a polygon's boundary.
y: the $y$ coordinates of a polygon's boundary.
ini.val: the list of initial values for the model parameters.
control: the list of control parameters for using the optim function in package stats.
par.list: the option of showing the list of parameters on the screen.
stand.fig: the option of drawing the observed and predicted polygons at the standard state (i.e., the polar point is located at (0, 0), and the major axis overlaps with the $x$-axis).
angle: the angle between the major axis and the $x$-axis, which can be defined by the user.
fig.opt: an optional argument of drawing the observed and predicted polygons at arbitrary angle between the major axis and the $x$-axis.
np: the number of data points on the predicted superellipse curve.
xlim: the range of the $x$-axis over which to plot the superellipse curve.
ylim: the range of the $y$-axis over which to plot the superellipse curve.
unit: the unit of the $x$-axis and the $y$-axis when showing the superellipse curve.
main: the main title of the figure.

Author

Peijian Shi pjshi@njfu.edu.cn, Johan Gielis johan.gielis@uantwerpen.be, Brady K. Quinn Brady.Quinn@dfo-mpo.gc.ca.

Details

The superellipse equation has its mathematical expression: $$r\left(\varphi\right) = a\left(\left|\mathrm{cos}\left(\varphi\right)\right|^{n}+ \left|\frac{1}{k}\,\mathrm{sin}\left(\varphi\right)\right|^{n}\right)^{-\frac{1}{n}},$$ where $r$ represents the polar radius at the polar angle $\varphi$. $a$ represents the semi-major axis of the superellipse; $k$ is the ratio of the semi-minor axis to the semi-major axis of the superellipse; and $n$ determines the curvature of the superellipse curve. This function is basically equal to the fitGE function with the arguments m = 4 and simpver = 9. However, this function can make the estimated value of the parameter $k$ be smaller than or equal to 1. Apart from the above three model parameters, there are three additional location parameters, i.e., the planar coordinates of the polar point ($x_{0}$ and $y_{0}$), and the angle between the major axis of the superellipse and the $x$-axis ($\theta$). The input order of ini.val is $x_{0}$, $y_{0}$, $\theta$, $a$, $k$, and $n$. The fitted parameters will be shown after running this function in the same order. The Nelder-Mead algorithm (Nelder and Mead, 1965) is used to carry out the optimization of minimizing the residual sum of squares (RSS) between the observed and predicted radii. The optim function in package stats was used to carry out the Nelder-Mead algorithm. When angle = NULL, the observed polygon will be shown at its initial angle in the scanned image; when angle is a numerical value (e.g., $\pi/4$) defined by the user, it indicates that the major axis is rotated by the amount ($\pi/4$) counterclockwise from the $x$-axis.

References

Shi, P., Gielis, J., Quinn, B.K., Niklas, K.J., Ratkowsky, D.A., Schrader, J., Ruan, H., Wang, L., Niinemets, Ü. (2022) 'biogeom': An R package for simulating and fitting natural shapes. Annals of the New York Academy of Sciences 1516, 123$-$134. tools:::Rd_expr_doi("10.1111/nyas.14862")

Shi, P., Huang, J., Hui, C., Grissino-Mayer, H.D., Tardif, J., Zhai, L., Wang, F., Li, B. (2015) Capturing spiral radial growth of conifers using the superellipse to model tree-ring geometric shape. Frontiers in Plant Science 6, 856. tools:::Rd_expr_doi("10.3389/fpls.2015.00856")

Examples

Run this code

data(whitespruce)

uni.C <- sort( unique(whitespruce$Code) )  
Data  <- whitespruce[whitespruce$Code==uni.C[12], ]
x0    <- Data$x
y0    <- Data$y
Res1  <- adjdata(x0, y0, ub.np=200, len.pro=1/20)
x1    <- Res1$x
y1    <- Res1$y

plot( x1, y1, asp=1, cex.lab=1.5, cex.axis=1.5, type="l", 
      col="grey73", lwd=2,
      xlab=expression(italic("x")), ylab=expression(italic("y")) )

# \donttest{
  x0.ini    <- mean( x1 )
  y0.ini    <- mean( y1 )
  theta.ini <- c(0, pi/4, pi/2)
  a.ini     <- mean(c(max(x1)-min(x1), max(y1)-min(y1)))/2
  k.ini     <- 1
  n.ini     <- c(1.5, 2, 2.5)
  ini.val   <- list( x0.ini, y0.ini, theta.ini, a.ini, k.ini, n.ini )
  Res2      <- fitSuper(x=x1, y=y1, ini.val=ini.val, unit="cm", par.list=FALSE,
                        stand.fig=FALSE, angle=NULL, fig.opt=TRUE,  
                        control=list(reltol=1e-20, maxit=20000), np=2000)
  Res2$par
  Res2$r.sq
# }

graphics.off()

Run the code above in your browser using DataLab