R/mod_QUAD.R

In MASSTIMATE: Body Mass Estimation Equations for Vertebrates

mod_QUAD <- function(xvar = extants$log.HCFC, yvar = extants$log.BM, 
                     xlab = "Log Stylopodial Circumference",ylab = "Log Body Mass",
                     xlim = c(1, 3.5), ylim = c(1.5, 8), cex = 0.5, pch = 4, col = "gray", 
                     fit.col = "black", pred.seq = seq(from = 1, to = 4, by = 0.01), 
                     PI.col = "blue", PI.lty = 2, CI.col = "red", CI.lty = 1, verbose = FALSE, plot = TRUE, ...) {
  extants <- data.frame(log.HCFC = c(2.795880017,2.911317442,1.342422681,1.366422957,1.662757832,1.679881942,1.911690159,2.021189299,1.937016107,1.707570176,1.70969387,1.651762447,1.725911632,1.511883361,2.113943352,1.767155866,1.916453949,1.900913068,1.886490725,1.873611197,1.935003151,2.352182518,2.171726454,2.089905111,2.333447274,2.26599637,2.102090526,2.10720997,1.482873584,1.45484486,1.361727836,1.382017043,1.700270937,1.630936119,1.42894429,1.788168371,1.873611197,1.703291378,1.669781615,1.663700925,1.167317335,1.2955671,1.401400541,1.405687787,1.467608106,1.783546282,2.122870923,2.142545884,1.808210973,2.29666519,2.435366507,2.46686762,1.656098202,1.503109437,1.359835482,1.245512668,1.407390904,1.880813592,1.788875116,2.049218023,1.938269483,1.365487985,1.544068044,2.270678836,1.777789187,1.55448916,1.679427897,1.213517757,1.416640507,1.340444115,1.381115081,1.508529719,1.48784512,1.92556991,1.456366033,1.805840549,1.373831145,1.147676324,1.178976947,1.056904851,1.152288344,1.297760511,1.01911629,1.008600172,1.077367905,1.28780173,1.290034611,1.146128036,1.227886705,1.122215878,1.098643726,1.252853031,1.224014811,1.685741739,1.372912003,1.393575203,1.353146546,1.826398782,1.28780173,1.447158031,1.322219295,1.095169351,1.1775365,1.152288344,1.170261715,2.117271296,1.533390708,1.716003344,1.478566496,1.767155866,1.69019608,1.550228353,1.521138084,1.024690862,1.191730393,1.117271296,1.096910013,1.091666958,1.08278537,1.486430479,1.400537989,1.764922985,1.457124626,1.217483944,1.431363764,1.286680969,1.728353782,1.664171705,1.442479769,1.545307116,1.277609214,1.378397901,1.204119983,1.322219295,1.267171728,1.051152522,1.062581984,1.2999429,1.200029267,1.737987326,1.434568904,1.553883027,1.722222464,1.089905111,1.315970345,1.136720567,1.376576957,1.536558443,2.01911629,1.946943271,2.104657791,1.202760687,2.127104798,2.123851641,1.977494969,2.555094449,2.401400541,2.492760389,2.068000226,2.247973266,2.028367884,2.236789099,2.33243846,2.076276255,2.273001272,2.378397901,2.18610838,2.214843848,1.886490725,2.292256071,2.071882007,2.439332694,2.555698895,2.290034611,2.461648568,2.314920056,2.155336037,2.351216345,1.889861721,2.173186268,2.245512668,2.407390904,2.411788005,2.439332694,2.608526034,2.413299764,2.318063335,2.56937391,2.683047038,2.622214023,2.190331698,2.350248018,2.302114377,1.717254313,2.16657811,1.921686475,1.827692289,1.808885867,2.241795431,1.904715545,1.258876629,1.301029996,1.268343914,1.478566496,1.2955671,1.131939295,1.258876629,1.170261715,1.231724383,1.128722284,1.049218023,1.155336037,1.324282455,1.633468456,1.324282455,1.200029267,1.516535374,1.549003262,1.200029267,1.719331287,1.576916956,2.11058971,1.804480189,1.706290957,1.753965866,1.301029996,1.135132651,1.396199347,1.491361694,1.408239965,1.27989498,1.367355921,1.446381812,1.135132651,1.139879086,1.152288344,1.176091259,1.393575203,1.539703239,1.64246452,1.529558673,1.482158695,1.66133934,1.328379603,1.452553063),
                        log.BM = c(6.548266545,6.808548551,2.439332694,2.648360011,3.562292864,3.56937391,4.133794301,4.505149978,4.187520721,3.544068044,3.681241237,3.361727836,3.77815125,3,4.781036939,3.924279286,4.1430148,4.248218561,3.986771734,3.944482672,4.129625452,5.307496038,4.851258349,4.785329835,5.361727836,5.161368002,4.63447727,4.672097858,3.075546961,3.041392685,2.653212514,2.584331224,3.677606953,3.385606274,2.855519156,3.770852012,4.218955418,3.806179974,3.789580712,3.812913357,2.161368002,2.414973348,2.84135947,3.041392685,3.054995862,3.995635195,4.633468456,4.62324929,3.84509804,5.309630167,5.638988159,5.650982244,3.414973348,3.021189299,2.431363764,2.250420002,2.625312451,4.058805487,3.544068044,4.456001437,4.167317335,2.414973348,2.977723605,5.190331698,3.806179974,3.076276255,3.361727836,2.245512668,2.556302501,2.618048097,2.785329835,2.855519156,2.979092901,4.096910013,3.079181246,3.951823035,2.664641976,2.012837225,2.036628895,1.977723605,1.806179974,2.29666519,1.826074803,1.84509804,1.785329835,2.371067862,2.274388796,2.025305865,2.008600172,1.934498451,1.892094603,2.041392685,2.212187604,3.348304863,2.574031268,2.544068044,2.602059991,3.956696565,2.224274014,2.69019608,2.322219295,1.826074803,1.973127854,1.908485019,1.913813852,4.752048448,3.111934276,3.653212514,3.18610838,3.531734309,3.633468456,3.278753601,3.024074987,1.748188027,2.079181246,1.897627091,1.792391689,1.84509804,1.785329835,3.079181246,2.666517981,3.910090546,2.740362689,2.376576957,3.053078443,2.278753601,3.795880017,3.478854968,2.883661435,3.044931546,2.537819095,2.67669361,2.201397124,2.537819095,2.290034611,1.723455672,1.781755375,2.414973348,2.123851641,3.92890769,2.84509804,3.176091259,3.875061263,1.770852012,2.371067862,1.886490725,2.77815125,3.049218023,4.260071388,4.356599436,4.579783597,2.028164419,4.781755375,4.691965103,4.342422681,6.071513805,5.491361694,5.736396502,4.469822016,5.067442843,4.405175546,5.035429738,5.409087369,4.492760389,5.071513805,5.498310554,4.88422877,4.881954971,4.14176323,5.128722284,4.706717782,5.478566496,5.832764305,5.092720645,5.718501689,5.352182518,4.914871818,5.200850498,4.385427515,4.820201459,5.120244796,5.435207103,5.544068044,5.418301291,5.995635195,5.491361694,5.397940009,6.101747074,6.278753601,6.156851901,4.956648579,5.397940009,5.143951116,3.488127496,4.469822016,3.792391689,3.653212514,3.737987326,5.224869595,4.343408594,2.413299764,2.643452676,2.7084209,2.785329835,2.709269961,2.298853076,2.411619706,1.995635195,2.298853076,2.049218023,1.982271233,2.356025857,2.608526034,3.431363764,2.572871602,2.385606274,3.08170727,3.13289977,2.123851641,3.704407927,3.041392685,4.857332496,3.936764998,3.684755622,3.730136004,2.729974286,1.897627091,2.744292983,3.200850498,2.88422877,2.267171728,2.662757832,3.252853031,2.08278537,2.139879086,2.238046103,2.285557309,2.798650645,3.072249898,3.451018452,3.171726454,2.983175072,3.583878598,2.444044796,3.139957756))
  fit <- lm(yvar~xvar+I(xvar^2)) #OLS regression analysis
  #for intervals
  seq <- pred.seq
  #prediction and errors
  quad.pred <- predict(fit, newdata = data.frame(xvar = seq), interval = "prediction")
  quad.conf <- predict(fit, newdata = data.frame(xvar = seq), interval = "confidence")
  if(plot) {
    #plot
    plot(xvar, yvar, xlim = xlim, ylim = ylim, xlab = xlab, ylab = ylab, cex = cex, pch = pch, col = col, ...)
    lines(seq, quad.pred[,1])
    #prediction intervals
    lines(seq, quad.pred[,3], col = PI.col, lty = PI.lty)
    lines(seq, quad.pred[,2], col = PI.col, lty = PI.lty)
    #confidence intervals
    lines(seq, quad.conf[,3], col = CI.col, lty = CI.lty)
    lines(seq, quad.conf[,2], col = CI.col, lty = CI.lty) 
  }
  if(verbose) {
    return(list(model = fit,
                prediction = quad.pred,
                confidence = quad.conf)) 
  }
}