Tree growth analysis issues, again, part five

Continuing in this series, the point of which is to explore issues in tree growth analysis and their relationship to claims made by Stephenson et al. (2014)

Previously, I derived a function relating tree above-ground biomass (AGB) to time, given an apriori equation from the literature (Chave et al., 2005) relating biomass to diameter that is similar to the one used by Stephenson et al. [I used this similar equation, rather than the one actually used by Stephenson et al. (2014), because doing so allowed for an easier derivation of the function relating diameter to time that would produce a constant mass growth rate.] The point of that exercise was to establish a reference frame of constant mass growth rate, from which I could then investigate what increasing and decreasing mass growth rates would imply for radial growth rates. The analysis showed me that models of that basic structure could only produce continuously increasing or decreasing mass (and radial) growth rates.

But continuously decreasing rates are biologically unreasonable, because they require impossible growth rates when trees are very young, and although continuously increasing rates do support what Stephenson et al. claim, they too are untenable. This is simply because, of course, nothing in nature can grow at an ever-accelerating rate; there has to be a deceleration at some point, at the very least to an asymptote, but more likely to an actual decline. The universe is not over-run by tree mass after all. In fact, one of the things that worries me about their paper is that their data shows very little rate deceleration at all, including in species that reach very large size. So…given the rate accelerations at the maximum tree sizes they report–the main point of their paper–then where are all the trees that are > 1.0X the mass of the largest they report? Shouldn’t they also exist? Did all their tree sizes just happen to fall short of the size threshold where growth rates start to decline? Or do they imagine that tree growth rates never slow down and maximum size is limited by various other forces? If so, what are they? Some things just don’t add up here.

With the previous model, all ring widths necessarily followed some type of exponential function with time (negative exponential whenever q >= 1.0, and positive exponential otherwise). But other patterns are quite possible, and in fact are routinely dealt with in tree ring studies, for example simple linear functions of negative slope. Another would be a unimodal response: rapid ring width increases in young trees, reaching a maximum sometime during the first few decades, before following an approximately negative exponential type of decline toward an asymptotic minimum thereafter. There are at least two functional types that can model this pattern with <= 3 parameters. These include Ricker, and “Holling Type IV” models. [Note that in some dendro studies, Rickers are referred to as “Huggershoff” functions].

So now, instead of working backward from mass as a function of diameter as before, I start with a (reasonable) ring width model and work forward. By “reasonable”, I mean one that displays a pattern that is commonly encountered in trees, whenever rings are observable. Because this ring width model is less restrictive than before, I can explore more relationships between growth rates and time. [I can also then use the actual biomass model that Stephenson et al used, instead of the simplified one, because I don’t have to solve anything, but still constrained by the biomass equation of Chave et al (2005) used by the authors, that relates mass to diameter.]

I’ll start with a Ricker function, which is a simple and flexible modification of a standard negative exponential function that allows asymptotes of zero (unlike e.g. the Holling Type 4), and having the form

y = c + ate^{-bt}

where y is ring width in mm, t is time in years, e is the base of natural logarithms (~2.72), and a, b, and c are model parameters. That’s a pretty nifty little equation if you look at it in terms of flexibility of response relative to structural complexity.

Here’s how I proceed. For sake of consistency with my previous analyses, I retain the mean radial (diameter) growth rates, each tree again reaching 120 cm dbh (mass = 20,635 kg) at t = 200 years. Within that constraint, I vary parameters a, b and c, which respectively determine the maximum ring width, the rate of change of ring width with time, and the minimum width. First, I choose the year in which maximum ring width occurs (t.max); I’ll test t.max = 33, 67 and 100. Parameter b, the “rate parameter”, is correlated with t.max, earlier t.max values having higher rates. Setting the first derivative to zero gives the relationship between t.max and b, which is simply b = 1/t.max, and so therefore b = 0.030, 0.015, and 0.010 respectively.

But I can also decrease (but not increase) the rate of ring width change with time via the scaling parameter, a, and the minimum value parameter, c. This I do as follows. The integral of the Ricker function from t = 1 to 200 years gives the radius, twice which is the diameter, which at year t = 200 must equal 120 cm. The diameter returned when a and c are non-influential (i.e., a = 1.0 and c = 0), for b = .05 at t = 200, is ~ 79.9 cm. Therefore parameter a needs to be scaled by a factor of 79.9/120, or 0.67. If we want to increase the minimum possible ring width, above c = 0, we can lower a accordingly, which will reduce the rate of ring width change. For example, for c = 1.0 mm minimum ring width, with a and b unchanged, then 200(1.0)/10 = 20 additional cm of growth, so a = (79.9 – 20)/120 = 0.50. We cannot go the other way though (reducing c while increasing a), because ring widths cannot be < 0 (trees don't shrink). I'll test values of c = 0.0 and 1.0 mm here, which when evaluated for each value of b mentioned above, gives 2 x 3 = 6 tested models. [Parameter “a”, being determined after the fact, does not influence the number of models tested.]

First, here’s the graph of the ring widths with time:
Fig1

which gives these diameters, and masses, as functions of time (notice that all lines meet at t=200 years, i.e. equal mean growth rates at that point):
Fig2
Fig3

These in turn give the following mass growth rates, first as a function of time, and second, as function of mass:
Fig4
Fig7

So, in these last two figures, we have results which are very different from what Stephenson et al. claim as their principal finding. Specifically, whether measured with respect to mass, or time, mass growth rates do not in any of the six models tested, increase continuously. Rather, they all hit a maximum rate at some intermediate time point, and then begin to decline. In some cases, this occurs early on and is very abrupt, while in others it occurs later and more gradually, but it always occurs, and no later than year 250.

Note that I’ve grown these trees over 300 years, 100 years beyond the point at which the mean growth rates are equal across all tested models. At year 200, two of the models are still showing increasing growth rates, because they’ve not yet hit their maxima. But by year 300 they have. Expanding, if I go back far enough in time, all of the models will give continuously increasing growth rates. Therefore, clearly, the length of time over which the measurements are made makes a difference in the conclusions to be drawn. And I could change the model parameters, beyond just the specific values I’ve tested here, but that result will still hold.

Someone might object at this point that the Ricker function is not the optimal one for modeling ring widths with time. For example, ring width patterns found in dendroclimatic studies often follow a simple negatively sloping linear function, or a mild negative exponential. While this is certainly true, the unimodal patterns modeled by the Ricker also found (hence the use of “Huggershoff” functions), but more importantly, a maximum growth rate followed by a decline, is a pattern not restricted to, or an artifact of, the Ricker function. This I will demonstrate in the next episode of our thrilling series on tree growth analysis. I’m told this is to be featured on Discovery Channel’s new series Real Men of Tree Growth Analysis.

Here be the R code once again:

## Mass growth rate forward modeling from Ricker-based ring width models  ##

density = 0.6						# mean AGB density, g/cm^3. Valued at 0.60 g/cm^3; Chave supplement gives values ranging from 0.54 to 0.8 roughly, but most are ~ 0.60
t.ref = 200; d.ref = 120				# reference point at which equal mean growth rates is true
time = 1:(1.5*t.ref)					# take to t = 300 years to see diverging dynamics

# Ricker function model y = ate^(-bt), with y = ring width and t = time
t.max = c(100/3, 2*100/3, 3*100/3)							# year at which max ring width occurs
b=1/t.max;b
c = c(0,1); comb = expand.grid(b,c)							# all combinations of b and c, given choices of t.max; full factorial design
colnames(comb)=c(&amp;quot;b&amp;quot;,&amp;quot;c&amp;quot;)
comb$add=2*(200*comb[,2])/10								# the diameter growth that is due strictly to parameter c whenever c != 0

ricker = function(a1=1, b1, c1=0,t=1:t.ref) c1 + a1*t*exp(-b1*t)		# Ricker function; 1st computation, a and c non influential
for(i in 1:nrow(comb)) comb$raw[i] = 0.2*sum(ricker(b1=comb[i,1]))	# What the resulting diameter would be (in cm)
comb$a = (120-comb$add)/comb$raw							# parameter &amp;quot;a&amp;quot; is determined such that a 200 year old tree is 120 cm dbh (and hence 20635 kg), the same as produced in the constant mass growth exercise (below)
comb$base = comb$raw * comb$a; comb							# what the resulting mass would be, in kg
form = rep(&amp;quot;ricker&amp;quot;,nrow(comb))
pars = data.frame(form,comb[,c(5,1,2)])						# final parameter table
rownames(pars) = paste(&amp;quot;model&amp;quot;,1:nrow(pars),sep=&amp;quot;&amp;quot;);pars
n.mods = nrow(pars)

ricker = function(a1,b1,c1,t=time) c1 + a1*t*exp(-b1*t)		# Ricker function; after parameters are set
results = array(dim=c(max(time),nrow(pars),7)); dim(results)
dimnames(results) = list(time, rownames(pars),
 c(&amp;quot;ring width (mm)&amp;quot;, &amp;quot;ln(dbh) (cm)&amp;quot;, &amp;quot;dbh (cm)&amp;quot;,
 &amp;quot;ln(mass) (kg)&amp;quot;, &amp;quot;mass (kg)&amp;quot;, &amp;quot;ln(mass growth rate) (kg/yr)&amp;quot;,
 &amp;quot;mass growth rate (kg/yr)&amp;quot;))
labels =  unlist(dimnames(results)[3])

# Use Chave's (2005) loglinear regression parameters relating ln(tree mass) to ln(diameter) for tropical moist forests:
int = -1.499; xcoef1 =  2.148; xcoef2 = 0.207; xcoef3 = -0.0281			# Note: adding xcoefs 2 and 3 give only very tiny improvements, based on Chave's AIC data in their Table 4.
for (i in 1:nrow(pars)) results[,i,1] = ricker(a1=pars[i,2], b1=pars[i,3], c1=pars[i,4])
results[,,3] = apply(2*results[,,1]/10, 2, cumsum)
results[,,2] = log(results[,,3]) 
results[,,4] = int + xcoef1*results[,,2] + xcoef2*results[,,2]^2 + xcoef3*results[,,2]^3 + log(density)
results[,,5] = exp(results[,,4])
results[-1,,6] = apply(results[,,4], 2, diff); results[1,,6] = NA 
results[-1,,7] = apply(results[,,5], 2, diff); results[1,,7] = results[1,,5]

# Plotting
set.seed(1); palette2 = sample(colors(),nrow(pars))
lgnd.text = paste(&amp;quot;a=&amp;quot;, round(pars[,2],2), &amp;quot;, b=&amp;quot;, pars[,3], &amp;quot;, c=&amp;quot;, pars[,4],sep=&amp;quot;&amp;quot;);lgnd.text
place = c(&amp;quot;topright&amp;quot;,&amp;quot;bottomright&amp;quot;,&amp;quot;topleft&amp;quot;,&amp;quot;bottom&amp;quot;,&amp;quot;bottom&amp;quot;,&amp;quot;topleft&amp;quot;,&amp;quot;bottomright&amp;quot;)

count = 0
for (i in c(1,3,5,7)){
 count = count +1
  windows()
  plot(time, results[,1,i], type=&amp;quot;l&amp;quot;, ylim=c(0,max(results[,,i])),
  main = paste(labels[i], &amp;quot;as function of time&amp;quot;), ylab=labels[i])
  for (j in 1:nrow(pars)) lines(time, results[,j,i], col = palette2[j], lwd=4)
  legend(place[count], legend=lgnd.text, fill = plotcolors, cex=1.5)
}

for (i in 3:5){
 count = count  + 1
 windows()
 plot(results[,1,i], results[,1,7], type=&amp;quot;l&amp;quot;,
  ylim = c(min(results[,,7]), max(results[,,7])),
  xlim = c(min(results[,,i]), max(results[,,i])),
  main = paste(&amp;quot;Mass growth rate as function of&amp;quot;, dimnames(results)[[3]][i], sep=&amp;quot; &amp;quot;),
  xlab = labels[i], ylab= labels[7])
  for (j in 1:nrow(pars)) lines(results[,j,i], results[,j,7], col = palette2[j], lwd=4)
  legend(place[count], legend=lgnd.text, fill = plotcolors, cex=1.5)
}
Advertisements

13 thoughts on “Tree growth analysis issues, again, part five

  1. Wow Jim…

    Haven’t been by in a couple weeks and look at you… Stephenson et al. has really got you all fired up. And to be quite clear – while I did get a copy of Stephenson (and thanks for pointing to the Open Access link!) I’ve not gotten through all of it. And I haven’t even gotten through all of your remarks yet. And even if I ever get that far I’m still likely to be only a touch more knowledgeable about these perennials (I’m an annual sort).

    So a couple questions:
    Is there some finite minimal size for the cambium of a tree? It would seem to me there must be. And if there is, and if every new growth season a new cambial layer is formed, then at some point a very old tree would seem to be faced with a physical requirement to produce an increasing amount of cambium every year just to make it around the old stem. So – the width of individual rings can decrease with time, but there comes a time where the overall DBH is sufficient that in order to get ’round the old you must make more total new (greater mass) than the previous year. If you can’t accomplish this for whatever reason then you are in a pickle – no? By this simplistic reasoning then, Stephenson et al. makes some sense – at least until a tree can no longer manage radial growth.

    Language has always fascinated me. So has anyone else noticed the term (perhaps in frequent use in forest ecology) used by Stephenson et al on the first page of their paper (second paragraph following abstract): log(tree mass) ?? I know, it’s math… but isn’t a significant part of a felled tree’s mass usually referred to as a log?

    Not a question but a ‘like’ – I like your including the R code.

    • Yes, I get fairly fired up when I see such grand claims in a glamour journal, especially when they’re based on questionable analyses, and use data transformations and figures that visually exaggerate already questionable findings (haven’t gotten to that latter point yet). There’s an enormous need for post-publication review in science, and it’s not even remotely close to getting done, except on blogs by people who very often come with a certain viewpoint rather than pure objectivity, in climate change and economics topics for example.

      I’m not 100% sure I’m following you, but the “minimum size” for the cambial layer is essentially a single cell layer. New cambial layers are formed each year only in the sense that this year’s cambial layer is not the exact same set of cells as last year’s–but it’s still in the same location, relative to the secondary xylem and phloem. As the tree gets bigger, that layer must contain more total cells to span the circumference, as you appear to state, which it does by a few cell divisions perpendicular to the radial divisions that form the wood rings (xylem) and phloem. Stephenson’s arguments are not affected by this issue though, because the cambium is an inconsequential fraction of the tree mass, and those few additional cells required to increase its circumference each year tinier still, and at any rate, their argument is about rates of mass growth relative to the mass of the tree. That is determined by the radial growth rate and the equation relating diameter to mass.

      Yes, a log is a section of the tree bole, typically a timber term. I also thought the unintentional play on words in their use of “log” transformations was kind of funny in that respect, and even moreso when you consider that the biomass equations they used for some of their North American species in fact only predicted the biomass of the tree bole, not the whole above ground tree.

  2. Isn’t a tree ring the result of each year’s new cambium? Maybe I missed something back in school. I was thinking each year’s cambial layer ended up as one more ring and from this all the ability to look at ring’s (ring width) as a proxy for environmental phenomena.

    So if ring width does have some biological relevance, how can it go to zero? This is what I was thinking of earlier – that there has to be some minimum greater than zero. In your model you have a couple of cases you test where you set C to zero (and if I read it correctly, C is the ring width right?). I get that ring width will tend to get narrower as the total stem gets larger – especially if there is no increase in the rate of growth. But it seems to me that it is at least mathematically possible for ring width to decline and yet the overall plant still realize a greater increase in total mass compared to the previous year. (the rate basis is one year, no?)

    I see in your reply that their basis is the mass of the tree. I had the impression it was the previous year’s incremental increase in mass… this is NOT a trivial difference and I need to do some checking here.

    • The ring is new, the cambium is not. The cambium is just a single layer of cells that regenerates itself constantly without differentiating into the xylem cells that make up the ring, or phloem cells (inner bark) to the outside of it. The rings are the proxies, never the cambium.

      Ring width can asymptote towards zero in environmentally limiting conditions (low temp, low light, drought), but rarely is it actually zero. Parameter c sets the minimum possible ring width. By itself it’s not that important–what’s important is the variation of the rings over time. Assuming generally higher growth rates, I could set the minimum width higher, but if the variation over time is similar to what I’ve presented here, I’m going to get similar growth rate results.

      Yes it is in fact possible to get decreasing ring widths and increasing mass growth rates as you state. I’ll demonstrate that more clearly in the next post, but you can see it even in this one: look at the end of the line for the model in which a=.07. Growth rate is increasing slightly while ring width is decreasing slightly.

      The rate basis is the mass of the tree at the start of the measured growth interval, which they state was typically a 5-10 year period for the tropical trees in their study. Which itself is a bit of a curiousity: time is the more logical and relevant basis. At any rate, they give similar, but not identical results.

  3. I’m no expert other than having spent decades with trees in the Adirondacks (upstate NY), California and now Colorado.

    I think you’re on to something here, Jim. I am simultaneously curious whether what you’re describing could relate to an effect that I’ve seen over time:

    In essence, many many trees are (for whatever reason) “unable” to retain 100% whole bark as they get older. Sometimes that’s due to external factors, but this post has me wondering…

    Could it be that physical limits to growth can easily lead to a tree more or less “splitting its seams” just like clothes that are too small? If it would require too much growth / too much mass to generate the next ring layer, wouldn’t it make sense for the outer bark to split?

    Perhaps some of the strip bark effects observed in BCP’s etc may begin with such anomalies.

    • The outer bark–the part that you see–definitely does split and crack, in the vast majority of species. That’s why the bark on the older parts of a tree very often looks rough and cracked whereas on younger stems and branches it is very often smooth. That results from the fact that the cambial layer that produces the outer bark, the cork cambium, is necessarily located a bit inside the surface, so the bark that’s more than a few years old gets pushed out to the surface, where the circumference is greater that it is at the cork cambium. So something has to give.

      With bristlecones and other high elevation species, strip bark is due largely to a different cause: abrasion. Basically the +/- constant direction of the prevailing winds scours the facing side of the tree, especially with winter ice, eventually abrading the bark and vascular cambium away, right down to the wood.

    • While much of strip bark effects may be due to abrasion, my observations on Almagre suggest other alternatives as well. In particular, it’s pretty obvious that severe storms can have a huge impact. The following photo makes this plain: https://picasaweb.google.com/Almagre.Bristlecones.2007/OtherInterestingTrees#5125159162286044882

      Have you seen studies of this abrasion effect? It may be testable with existing data (climateaudit.info/data/colorado/2007_Almagre_Info.xls — Bark column W in the Trees tab.)

      Below is a set of compass-measured angles of Strip or Whole bark on a good-sized set of BCP’s at Almagre, all one site within a short distance of one another. I’ll see if I can work out a way to visually graph this information, and possibly analyze it to see the extent of correlation to prevailing wind/ice-driven abrasion. Interesting!

      Tree ID Bark (Whole or Strip) coverage, degree ranges
      024 W0-360
      026 W0-360
      028 W10-200 @gnd
      030 W60-70, W100-130
      031 W130-260, SW side
      033 W40-250 slight spiral
      037 W30-140
      047 S250-280 nr base
      048 S200-330 spiral
      040 S260-360
      041 W~10-180
      042 W20-180
      044 W60-210
      045 W100-260
      032 W150-170
      049 W~20-200
      050 S0-360
      051 S~200-270
      052 W0-360
      053 S0-360
      054 S~180-190
      055 W~0-180; 3′ up, S320-0
      056 S~270-340
      057 W~90-240
      058 S0-360
      059 S~10-30, S~50-80, S~160-230, S~260-310
      060 S~0-10, S~220-300
      061 S0-360
      062 W120-160
      063 W40-170
      064 W100-140, W230-235
      065 W80-200
      066 W20-240
      068 W0-360

    • You could be right–I was just giving the common explanation, which probably holds in most cases. Will check the collection to see what I have.
      Interesting data, but I can’t understand how you can get something like S0-360 or W40-170, which seem to be internally contradictory given the definitions.

    • S0-360 is a dead tree. There’s no bark left.
      W40-170 says there’s bark from 40-170 degrees on the compass dial — so 130 degrees of coverage from ~NE to ~S side of the tree.

      That’s the entire set of bark data I had time to collect before the weather overwhelmed me. All of this is inspiring me to get back up there and do some additional measurements. (One thing I’m hoping for, in light of improved technology, is to create 3D images of trees/bark. Some new tools show possibility of doing that with lightweight, portable equipment. Soon possibly even with just a special smart phone!)

    • Anything that inspires one to get into a subalpine tree stand should be vigorously pursued!

      Here’s a study you’ll like that relates exactly.

      The confusing part about the labels is that W is supposed to represent whole bark and S is strip bark, so I still don’t see it. If there’s bark for 130 degrees of coverage, that’s not whole bark then right?

    • Oh! I’m sorry, that’s a good catch nobody else noticed. I’ll fix the notes. “W” means “bark present” and “S” means “bark stripped” (Sometimes it was much easier to measure what was left vs what was gone…)

    • Thanks for the link. Interesting study! I think I can replicate their presentation with “our” trees, at ten (instead of five) degree intervals.

      Their results make intuitive sense. I’ve lived among similar trees at 7000 ft for the last two decades. What they don’t mention (but is certainly compatible with their conclusion): wind-driven “particulate abrasion” can be anything from small sand particles, to branches or entire tree trunks in a big storm. My observation of the SB BCP’s is that much of the stripping appeared to be due to violent physical processes. As noted by Steve McIntyre back in 2007, there was a major 19th century storm in the CO Front Range; the fact that the storm correlates well with later SB growth surges is at least suggestive that storm-related damage played a significant role in stripping many trees.

      Well, that’s some good food for thought, and inspiration for a future adventure. Please give a yell if you ever want to visit some cool trees. I can see Almagre from our deck 🙂

      Blessings,
      Pete

    • Yeah, I thought they did a good job of assessing all the various possibilities, including sun scald, abrasion, storm damage, etc. Good paper.

      If I’m in the Front Range with some time to spare I’d love to go up there with you.

      Get yourself out to the Snake Range in E NV, Wheeler Peak in particular, if you get a chance and haven’t been yet. Spectacular setting and P. longaeva groves on that mountain. In fact the Prometheus tree (oldest tree ever known) is there and is a fun though depressing excursion (it was cut down in 1964).

Have at it

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s