forked from RemkoDuursma/duursma2015plosone
-
Notifications
You must be signed in to change notification settings - Fork 0
/
manuscript.Rmd
333 lines (208 loc) · 28.6 KB
/
manuscript.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
---
title: "plantecophys - an R package for analysing and modelling leaf gas exchange data"
author: ''
date: ''
output:
word_document:
reference_docx: manuscriptstyle.docx
pdf_document: default
html_document:
number_sections: yes
csl: plos-one.csl
bibliography: references.bib
---
\
Remko A. Duursma^1\*^
\
^1^ Hawkesbury Institute for the Environment, Western Sydney University, Penrith, NSW, Australia
\
\* *Corresponding author*
\
Email: [email protected]
\
*Short title*: The plantecophys R package
```{r echo=FALSE, message=FALSE, warning=FALSE}
# Load packages and custom functions
source("load.R")
# Find current SHA of plantecophys (gives '???' on fail)
cur_sha <- "???" #get_git_SHA("c:/repos/plantecophys")
# Analysis for Figure 5.
source("R/analysis.R")
# Make PDFs
source("R/make_figures.R")
# This document can be compiled with the 'Knit Word' button in Rstudio, or:
# library(rmarkdown)
# render("manuscript.Rmd", "word_document", "manuscript.docx")
# which assumes you have Pandoc installed (unless you run this from Rstudio),
# and that you have an internet connection.
```
# Abstract {.unnumbered}
Here I present the R package 'plantecophys', a toolkit to analyse and model leaf gas exchange data. Measurements of leaf photosynthesis and transpiration are routinely collected with portable gas exchange instruments, and analysed with a few key models. These models include the Farquhar-von Caemmerer-Berry (FvCB) model of leaf photosynthesis, the Ball-Berry models of stomatal conductance, and the coupled leaf gas exchange model which combines the supply and demand functions for CO~2~ in the leaf. The 'plantecophys' R package includes functions for fitting these models to measurements, as well as simulating from the fitted models to aid in interpreting experimental data. Here I describe the functionality and implementation of the new package, and give some examples of its use. I briefly describe functions for fitting the FvCB model of photosynthesis to measurements of photosynthesis-CO~2~ response curves ('A-C~i~ curves'), fitting Ball-Berry type models, modelling C3 photosynthesis with the coupled photosynthesis-stomatal conductance model, modelling C4 photosynthesis, numerical solution of optimal stomatal behaviour, and energy balance calculations using the Penman-Monteith equation. This open-source package makes technically challenging calculations easily accessible for many users and is freely available on CRAN.
\
Keywords: photosynthesis, stomatal conductance, portable gas exchange instrument, process-based model, software
# Introduction
Since the advent of portable gas exchange instruments [@mcdermitt1989; @leuning1989], a wealth of data on leaf gas exchange of CO~2~ and H~2~O has been collected [@long2003]. These data play a central role in physiological plant ecology [@pearcey1989], to better understand and quantify inter-specific differences in photosynthesis and transpiration, and to quantify and model the rapid response to changes in environmental drivers such as light, humidity and temperature. Not only do leaf gas exchange data allow detailed studies of the underlying plant physiology, they are also used to parameterize an important component of process-based models of vegetation function used to predict global water and carbon cycling [e.g. @bonan2014; @dekauwe2015].
\
The photosynthesis model of Farquhar, von Caemmerer and Berry [@farquhar1980] (the 'FvCB model') is widely used in interpreting and modelling leaf gas exchange, by providing comparable metrics of the photosynthetic capacity, and predicting the response of photosynthesis to changes in the CO~2~ concentration inside the leaf air space (C~i~). This widely cited model is embedded in many process-based models of vegetation function [@sellers1997; @bonan2014]. The key prediction of the model is the response of photosynthesis to [CO~2~] inside the leaf (either chloroplastic [CO~2~], C~c~, or intercellular [CO~2~], C~i~). It can also account for changes in leaf temperature if the various temperature sensitivities are parameterized [@caemmerer1981; @bernacchi2001; @medlyn2002]. To employ the model, it is generally fit to observations of net photosynthesis along a range of [CO~2~] concentrations, yielding well-known measures of photosynthetic capacity (V~cmax~ and J~max~, and optionally R~d~) [@wullschleger1993].
\
I do not repeat a detailed description of the FvCB model here, as it has been described many times [e.g. @bernacchi2001; @medlyn2002]. But generally it is of the form,
\
(1)
$$A_n = min(A_c, A_j) - R_d$$
\
where A~n~ is the net rate of CO~2~ assimulation, A~c~ is the gross photosynthesis rate when Rubisco activity is limiting, A~j~ when RuBP-regeneration is limiting, and R~d~ the rate of dark respiration (see Fig. 1a). A~c~ and A~j~ are non-linear functions of the chloroplastic CO~2~ concentration (C~c~), both of the form $k_1(C_c{\Gamma}^{∗})/(k_2 + C_c)$, where ${\Gamma}^{∗}$ is the CO~2~ compensation point without R~d~, and k~1~ and k~2~ are different parameter combinations for A~c~ and A~j~. The details of these functions and the temperature dependence of the various parameters are described elsewhere [@medlyn2002].
\
In the practical application of the FvCB model, when leaf gas exchange is measured with a portable gas exchange instrument, estimates of C~c~ are difficult to obtain because they require an estimate of the mesophyll conductance (g~m~). In this case, it is customary to us the intercellular [CO~2~] concentration (C~i~) as the driver of photosynthesis. This approach is useful because C~i~ can be estimated from concurrent measurements of CO~2~ and H~2~O flux [@sharkey1982]. In the remainder of this article I will use C~c~ as the driver of photosynthesis, but point out that this can be replaced by C~i~ if the user does not have an estimate of g~m~. When g~m~ is known, C~c~ is calculated from C~i~ with Eq. 2.
\
(2)
$$C_c = C_i - A_n/g_m$$
\
where g~m~ is the mesophyll conductance (mol m^-2^ s^-1^). Although this method assumes that g~m~ is constant for a given leaf, it is well known that g~m~ responds dynamically to fluctuations in environmental drivers [@flexas2008] although some of the variation in g~m~ may due to artefacts related to (photo-)respiratory effects on measured g~m~ with standard methods [@tholen2012]. Because no model has been developed to date that adequately captures the variation in g~m~, it is a constant parameter in the implementation presented here. However, it can still be used to study the effects of non-constant g~m~ on rates of photosynthesis and its response to environmental drivers, as the parameter can be varied in model simulations.
\
**Figure 1**. (a) Leaf photosynthesis - CO~2~ response curve as modelled with the FvCB model. (b) The intersection of the supply and demand curves of photosynthesis. The `Photosyn` function solves for C~i~ if g~s~, V~cmax~, J~max~ and R~d~ (and other parameters to the FvCB model) are known.
\
Fitting the FvCB model to data requires some finesse because net photosynthesis is modelled as a minimum function of two non-linear equations that is sometimes difficult to fit. Moreover, sample sizes collected are often small due to time constraints in the field. A widely used published method requires the user to specify the transition of V~cmax~ to J~max~ limitation [@sharkey2007], a process that is both arbitrary and prevents batch analysis. Another method [@leafweb] requires online submission of data and fits the model without much control or knowledge of the fitting process (following [@gu2010]), and does not report standard errors of the estimated parameters. Undoubtedly many more implementations of the fitting process have been developed over the years, but few of these are made publicly available (but see available online tools [@dekauwefitfarquhar; @landfluxtools]). What is missing is an open-source tool that can be used for reproducible and transparent analysis of A-C~i~ curves.
\
Through Eq. 1, we have a dependency of photosynthesis on the availability of the substrate, C~c~. To estimate C~c~ itself, we need C~i~, which can be estimated when with stomatal conductance to CO~2~. From Fick's law, we can relate A~n~ to g~s~ and C~i~ as,
\
(3)
$$ A_n = \frac{g_s}{1.6}(C_a - C_i) $$
\
where g~s~ is the conductance to H~2~O (the factor 1.6 converts to conductance to CO~2~). We now have two equations for A~n~: the 'demand function' (Eq. 1), and the 'supply function' (Eq. 3). At steady state these two equations should be equal, which can be graphically shown as in Fig. 1b (cf. [@farquhar1982]).
\
Because g~s~ itself responds to environmental drivers, another expression is needed to end up with a fully coupled model of leaf gas exchange. The most widely-used, though empirical, g~s~ model is the Ball-Berry [@ballberry1987] class of models. This model posits an entirely empirical equation that describes the response of g~s~ to air humidity, CO~2~ and A~n~. This way, effects of leaf temperature and PPFD - both of which are known to affect g~s~ - are modelled through the dependency of A~n~ on these drivers. A general form of the Ball-Berry model is,
\
(4)
$$ g_s = g_0 + g_1 \frac{A_n}{C_a}f(D)$$
\
where D is the vapour pressure deficit (kPa), g~0~ and g~1~ are empirical parameters, and f(D) can be one of many functions that describe the response to the vapour pressure deficit (D, [@leuning1995; @medlyn2011]) or relative humidity [@ballberry1987]. An alternative approach to modelling g~s~ is through the hypothesis that stomata act optimally in the sense that they maximize photosynthesis while minimizing water loss. This hypothesis was first developed by Cowan and Farquhar [@cowan1977] and has seen many applications. Medlyn et al. [@medlyn2011] showed that the optimality hypothesis, when coupled to the FvCB model, leads to an expression analogous to the Ball-Berry type models (Eq. 4), but with a different D response function (f(D) in Eq. 3) compared to the original Ball-Berry model.
\
Finally we can combine the biochemical demand function of photosynthesis (Eq. 1) with the supply function (Eq. 2) and an expression for the dependency of g~s~ on environmental drivers (Eq. 3). This 'coupled' leaf gas exchange model [@leuning1990; @collatz1991; @leuning1995] is implemented in many process-based ecosystem and global land surface models [@wang1998; @bonan2014; @duursma2012; @sellers1997]. This model allows prediction of A~n~, g~s~ and leaf transpiration rate in response to all major environmental drivers (except soil water limitation), and incorporates key leaf traits (g~1~, V~cmax~, J~max~, R~d~, and their temperature dependencies).
\
Despite the widespread use of the FvCB model and the coupled leaf gas exchange model, tools to analyse data and perform simulations are scattered and subject to little standardization. Fitting the FvCB model to CO~2~ response curves is a standard procedure but different methods can yield different parameter values, making comparisons difficult. The coupled leaf gas exchange model is not straightforward to implement, and I do not know of any standalone open-source implementations. I here describe the plantecophys package, implemented in the R language [@rfoundation]. The code is freely available (without restrictions), and managed with a version control system. The package is the result of our work on leaf and canopy modelling of photosynthesis and stomatal conductance [@medlyn2011; @barton2012; @peltoniemi2012; @medlyn2013; @duursma2013; @duursma2014; @lin2015], with many additions based on user requests.
# Design and implementation
## The main functions
The main tools included in the plantecophys package are to a) fit A-C~i~ curves to estimate V~cmax~, J~max~ and R~d~, b) fit Ball-Berry type models, c) simulate from the coupled leaf gas exchange model and d) calculate the optimal stomatal conductance. The key functions in the package are summarized in Table 1.
```{r results='asis', echo=FALSE, message=FALSE, warn=FALSE}
f <- c("`fitaci`","`fitBB`","`FARAO`","`Photosyn`","`PhotosynEB`","`AciC4`", "`RHtoVPD` etc.")
d <- c("Fit, summarize, plot and simulate photosynthesis-[CO~2~] response curves (A-Ci curves)",
"Fit Ball-Berry type models of stomatal conductance",
"Estimate optimal stomatal conductance with a numerical implementation of the Cowan-Farquhar hypothesis",
"Simulate C3 photosynthesis and transpiration with the coupled leaf gas exchange model. Also simulates the FvCB model when either C~i~ or g~s~ is given as input.",
"Estimate leaf temperature from energy balance, when a significant leaf boundary layer is present",
"Simulates the dependence of C4 photosynthesis on the intercellular CO~2~ concentration",
"Convert between commonly used units (relative humidity, vapour pressure deficit, dewpoint temperature)")
df <- data.frame(Function=f, Description=d)
suppressWarnings(print(ascii(df, caption="Table 1. Main functions in the plantecophys package.",
include.rownames=FALSE),
type='pandoc'))
```
## Language
The 'plantecophys' package is implemented in R, has no dependencies on other packages, and does not require compilation (i.e. it is written in native R only). As such it builds easily, and is highly portable. The source code is maintained with git version control, and is hosted in an online repository (http://www.bitbucket.org/remkoduursma/plantecophys), from which a development version of the package can easily be installed. The repository includes an issue tracker, where users can suggest changes or report bugs. This paper describes version `r packageVersion("plantecophys")` (git SHA `r cur_sha`).
\
All code used in this article (including the code to generate the article written in markdown, all figures and full example code), can be downloaded from http://www.github.com/duursma2015plosone. The repository also includes code to demonstrates how to extract additional statistics from fitted A-C~i~ curves.
# Results and Discussion
### Fitting A-C~i~ curves
The `fitaci` function fits the FvCB model, yielding estimates of V~cmax~, J~max~ and R~d~ and their standard errors. Instead of fitting the minimum function (Eq. 1), `fitaci` fits the hyperbolic minimum of A~c~ and A~j~, which avoids a discontinuity (Eq. 5).
\
(5)
$$A_m = \frac{A_c+A_j - \sqrt{(A_c+A_j)^2-4 \theta A_c A_j}}{2 \theta} - R_d$$
\
where $\theta$ is a shape parameter, set to 0.9999, and A~m~ is the hyperbolic minimum of A~c~ and A~j~. The fit of the FvCB model to data is achieved with non-linear least squares, and standard errors of the parameters are estimated with standard methods (`nls` function in base R, see [@ritz2008]). The `fitaci` function includes methods to estimate appropriate starting values from the data, and attempts the fits along a wide range of possible starting values. Optionally, R~d~ can be provided as a known value, otherwise it is estimated from the A-C~i~ curve. The user does not have to provide the transition point (see Fig. 1), as this is estimated by `fitaci` automatically. It is however an option to fix the transition point (via the `citransition` argument), which may be helpful to check whether the best fit was achieved. Finally, the user can provide an estimate of mesophyll conductance (g~m~) (following [@ethier2004]), in which case the fitted values of V~cmax~ and J~max~ can be interpreted as chloroplastic rates.
\
Because the fitting uses non-linear least squares, standard methods can be employed to estimate standard errors (SE), confidence intervals, and correlation of the fitted parameters. The `fitaci` function returns by default the SE and confidence intervals, and the built-in help page for the `fitaci` function shows how the `nlstools` package can be used to provide a detailed overview of the statistics of the non-linear least squares fit.
\
Required inputs are measurements of A~n~ and C~i~, and optionally leaf temperature (T~leaf~), and photosynthetically active radiation (PAR). Also required are estimates of Michaelis-Menten constants (K~c~, K~o~ or the combination K~m~) and $\Gamma^{*}$. In the FvCB model, J~max~, V~cmax~ and leaf respiration (R~d~) (and other parameters like $\Gamma^{*}$, K~c~ and K~o~) all depend non-linearly on T~leaf~. The `Photosyn` function incorporates standard temperature sensitivities for all parameters of the FvCB model (following [@medlyn2002]). Optionally, measured (or otherwise modelled) K~m~ and $\Gamma^{*}$ can be provided as input.
\
The function takes a dataframe as input, which includes measurements of A~n~, C~i~ and optionally T~leaf~ and PAR, and is easily used like this:
```
# Fit FvCB model
f <- fitaci(mydata)
# Print a summary with coefficients and more
f
# Make standard plot
plot(f)
```
\
The output of the above example is shown in Fig. 2. Additionally, the batch utility `fitacis` can be used to fit many curves at once, for example one for each species or site in a dataset. I show this functionality in the example application further below.
```{r eval=TRUE, echo=FALSE, message=FALSE}
# acidata1 is an example dataset in plantecophys (see ?acidata1)
acidata1$PPFD <- 1800
exfit1 <- suppressWarnings(fitaci(acidata1))
fitc <- exfit1$pars
```
\
**Figure 2**. Standard output from the `fitaci` function. A~n~ is the net photosynthetic rate, C~i~ the intercellular CO~2~ concentration. Symbols are measurements, the black line the fitted FvCB model of photosynthesis. Colored lines indicate the two photosynthesis rates in the FvCB model. In the default mode, the `fitaci` function estimates V~cmax~, J~max~ and R~d~ from the fitted curve. Optionally, R~d~ is provided as an input, for example when it was measured separately. In this example, V~cmax~ was estimated as `r round(fitc["Vcmax","Estimate"],1)` (SE `r round(fitc["Vcmax","Std. Error"],2)`), J~max~ was `r round(fitc["Jmax","Estimate"],1)` (SE `r round(fitc["Jmax","Std. Error"],2)`) and R~d~ was `r round(fitc["Rd","Estimate"],1)` (SE `r round(fitc["Rd","Std. Error"],2)`). Assumed parameters were K~m~ = `r round(exfit1$Km,0)` and ${\Gamma}^{∗}$ = `r round(exfit1$GammaStar,1)` (all in units of $\mu$mol m^-2^ s^-1^). The R^2^ of a regression of measured vs. fitted was 0.99.
\
A C4 model of leaf photosynthesis [@caemmerer2000] is also implemented (in `AciC4`), but at the moment it is only possible to fit the C3 model of leaf photosynthesis to A-C~i~ curves.
### Fitting stomatal conductance models
The straightforward `fitBB` function provides an interface to non-linear or linear regression to fit one of three stomatal conductance models [@ballberry1987; @leuning1995; @medlyn2011]. This yields estimates of g~1~ and (optionally) g~0~, which are necessary inputs to the coupled leaf gas exchange model. Note that the user must provide stomatal conductance to H~2~O (not CO~2~) as input to the fitting process, which is the standard output of portable gas exchange instruments. This function is demonstrated in the example application further below.
### Coupled leaf gas exchange model
The intersection of the supply and demand curves of photosynthesis (Fig. 1b) gives the steady-state intercellular CO~2~ concentration (C~i~). This is solved by the `Photosyn` function. This flexible interface can be used to either 1) estimate A~n~ when C~i~ is known (`Photosyn(Ci=...)`; equivalent to `Aci(...)`), 2) estimate A~n~ when g~s~ is known (`Photosyn(GS=...)`) (cf. Fig. 1b) or c) solve for C~i~ from the coupled leaf gas exchange model (Eqs. 1,3,4).
\
To demonstrate the use of the coupled gas exchange model, I visualize the temperature response of A~n~ when both T~leaf~ and D are varying. In field conditions, D is always strongly positively related to T~leaf~. The consequence is that when studying D or T~leaf~ responses in the field, both drivers have to be accounted for simultaneously [@lin2012; @duursma2014]. Figure 3a shows simulated A-C~i~ curves and the solutions of the coupled leaf gas exchange models at a range of T~leaf~ and corresponding D (calculated following [@duursma2014]). Both V~cmax~ and J~max~ have a peaked response to T~leaf~, so that at a given C~i~, A~n~ first increases with T~leaf~ and then decreases (lines, Fig. 3a). As a result of increasing D, the modelled C~i~ decreases (symbols, Fig. 3a, as a consequence of Eq. 4). The net result is a peaked response of A~n~ as a function of D (Fig. 3b).
\
**Figure 3**. Response of A~n~ and C~i~ to combined changes in T~leaf~ and D. (a) Lines are A-C~i~ curves simulated at a range of values for T~leaf~. Symbols are the solutions of the coupled leaf gas exchange model, while also taking into account the correlation between D and T~leaf~ (based on an empirical relationship [@duursma2014] : D = 0.000605*T~air~^2.39^). Note that as T~leaf~ and D increase, C~i~ decreases. (b) The corresponding temperature optimum of A~n~. Symbols are the same as in panel (a) but plotted against T~leaf~.
\
The simplified code to produce Figure 3b, using the `Photosyn` function, is given below. Note that in this example the default values of many parameters (e.g. J~max~, g~1~) are used in the call to `Photosyn`, but all of these can be set by the user.
```
# Set range of leaf temperature
tleafs <- seq(5, 40, by=5)
# Define D as a function of Tleaf
vpdfun <- function(tair)0.000605*tair^2.39
# Simulate.
run1 <- Photosyn(Tleaf = tleafs, VPD=vpdfun(tleafs))
# Plot (produces Figure 3b minus the special formatting)
with(run1, plot(Tleaf, ALEAF))
```
\
The `Photosyn` function assumes that the boundary layer conductance (g~bl~) is high compared to g~s~, so that T~leaf~ is close to T~air~. As an alternative, the `PhotosynEB` function calculates T~leaf~ from the leaf energy balance. Transpiration is calculated with the Penman-Monteith equation [@jones1992], which accounts for boundary layer effects. The details of `PhotosynEB` are not described here (see the built in help file for more information), because it is very similar to other implementations [@wang1998; @buckley2014].
### Numerical solution of optimal stomatal conductance
The `FARAO` function (FARquhar And Optimality) calculates optimal stomatal conductance based on the Cowan-Farquhar [@cowan1977] hypothesis that stomata respond to environmental drivers in order to maximize photosynthesis while minimizing water loss. This implementation was used by Medlyn et al. [@medlyn2011] to compare a simplified model of optimal stomatal conductance to the full numerical solution.
To find optimal stomatal conductance, `FARAO` finds the C~i~ for which the quantity $A_n - \lambda E$ is maximal, where E is the leaf transpiration rate and $\lambda$ is the marginal cost of water (an empirical parameter related to g~1~, see [@cowan1977; @medlyn2011]). A~n~ is calculated directly as a function of C~i~ via the FvCB model (Eq. 1), g~s~ is calculated by rearranging Eq. 3, and E is calculated assuming perfect coupling (thus $E = g_s D/P_a$, where P~a~ is atmospheric pressure). This numerical routine does not need specification of an f(D) function as in Eq. 4, instead, this function is an emergent property. In Fig. 4a, I have calculated $A_n - \lambda E$ across a range of C~i~ values, and for different values of D. The `FARAO` function calculates the optima of these curves, which can for example be used to study the response of stomatal conductance to D (Fig. 4b).
\
**Figure 4**. Visualization of the optimal model of stomatal conductance. Provided we have an estimate of the 'cost of water' ($\lambda$, mol C mol H~2~O^-1^), stomata act to maximize photosynthesis minus transpiration. In (a), individual curves at a range of values for the vapour pressure deficit (D) are plots of $A-\lambda E$ as a function of C~i~, demonstrating that an optimum C~i~ exists. The `FARAO` function finds this optimum numerically and calculates corresponding A~n~ and g~s~. The corresponding response of g~s~ to D is shown in panel (b).
\
Optionally, the `FARAO` function accounts for the presence of a leaf boundary layer (when `energybalance=TRUE`). In that case it uses `PhotosynEB` (see description above) to calculate A~n~ and E, and solves for T~leaf~. A very similar method was employed by Buckley et al. [@buckley2014], who demonstrated that when a boundary layer is present, frequently an optimal g~s~ cannot be found.
### An example application
To demonstrate a practical application of the key functions in the package, I use field-collected data from Medlyn et al. [@medlyn2007; @medlynfigshare] on *Eucalyptus delegatensis*. Both A-C~i~ curves and 'spot gas exchange' data (i.e. leaf gas exchange measurements at prevailing environmental conditions) were collected. Using the `fitacis` function, it is straightforward to fit all `r length(acifits)` curves to the A-C~i~ data, and make standard plots of the fitted curves (shown in Fig. 5a). The fitted coefficients can be extracted using the `coef` function, and used to plot a comparison of fitted V~cmax~ and J~max~ values, which show the typical correlation between the two (Fig. 5b).
\
Next, I fit Eq. 4 to the spot gas exchange data, yielding an estimate of g~1~ (Fig. 5c). In this example, I used the model of Medlyn et al. [@medlyn2011], which is given by Eq. 6 (in this example, I assumed g~0~ = 0)
\
(6)
$$ g_s = g_0 + 1.6(1 + \frac{g_1}{\sqrt(D)})\frac{A_n}{C_a}$$
In Fig. 5c, modelled g~s~ is compared to measurements (Fig. 5c). To compare the model prediction of instantaneous transpiration efficiency (A~n~/E) to measurements along the variation in D (Fig. 5d), Eq. 6 can be rearranged to give (cf. [@duursma2013], where it is assumed that g~0~ = 0)
\
(7)
$$A_n/E = \frac{C_a P_a}{1.6 (g_1 D_s^k + D_s)}$$
Because `fitBB` can fit a number of Ball-Berry type variants, the various models can be easily compared in terms of goodness of fit. This simple example application is available in the published repository (see Methods), and simplified code for this example (panels a-c only), omitting special formatting and a few minor settings, is given below.
```
# Fit A-Ci curves.
# In this case, each separate curve is indexed by a column named 'Curve',
# and the data were already read into a dataframe (tumh)
acifits <- fitacis(tumh, "Curve")
# Plot all A-Ci curves in one panel, highlight one fitted curve.
plot(acifits, "oneplot", highlight="25")
# Fit Medlyn et al.'s (2011) version of the Ball-Berry model
# Data are already read into a dataframe (tumspot),
# and have standard names (or they can be set).
gfit <- fitBB(tumspot, gsmodel="BBOpti")
# Plot measured versus modelled, by predicting from the fitted model.
tumspot$GSpred <- predict(gfit$fit, tumspot)
with(tumspot, plot(GSpred, Cond))
```
\
**Figure 5**. Example application of the plantecophys package to A-C~i~ curves and spot gas exchange measurements on *Eucalyptus delegatensis* (data can be downloaded from [@medlynfigshare]) (a) Fitted A-C~i~ curves with one curve highlighted (b) Estimates of J~max~ plotted against V~cmax~, obtained from the fitted curves in panel (a). Solid line is a regression line (J~max~ = `r round(lmjvt$estimate[1],2)` + `r signif(lmjvt$estimate[2],2)` V~cmax~, R^2^ = `r signif(lmjvg$r.squared,2)`) with a 95% confidence interval for the mean. (c) Modelled (with the model of Medlyn et al. 2011) versus measured g~s~ (p `r pval(lmgpredg$p.value)`, R^2^ = `r signif(lmgpredg$r.squared,2)`). Measurements included a wide range of environmental conditions (PAR, T~leaf~, D). In this example, only g~1~ was fit (estimate = `r signif(g1,3)`, 95% CI = `r paste(signif(g1ci,3), collapse=" - ")`). (d) The predicted response of ITE (A~n~ / E) as a function of D from the fitted model in panel (c) (solid line), and the measurements from panel (c) when PAR > 1000.
# Conclusions
We need an open source set of tools to analyse leaf gas exchange data, as these data form a cornerstone of plant physiological ecology. At the moment there are no publicly available tools to fit A-C~i~ curves or perform simulations with the coupled leaf gas exchange model that can be used as part of a reproducable workflow. The plantecophys R package is implemented in widely used language for data analysis. The package includes a useful set of tools to perform standard, and more advanced, analyses of leaf gas exchange data. The open source framework combined with version control allows further development of the code.
# Availability and requirements
**Project name**: plantecophys
**Project Stable Release**: cran.r-project.org/package=plantecophys
**Project Home Page**: http://www.bitbucket.org/remkoduursma/plantecophys
**Project Issue Tracker**: http://www.bitbucket.org/remkoduursma/plantecophys/issues
**Operating System(s)**: Platform Independent
**Programming Language(s)**: R
**Other Requirements**: none
# Acknowledgments {.unnumbered}
David Ellsworth, Belinda Medlyn and Martin De Kauwe are acknowledged for discussions on fitting A-C~i~ curves. Special thanks to Belinda Medlyn for the code in the MAESTRA model which provided the basis for an early version of the `Photosyn` function, and for sharing the data used in the example application. Thanks to John Drake for testing and suggesting new features.
# References {.unnumbered}