Various aspects of uncertainty have become topical in pest risk modelling discussions. A recent contribution to the literature sought to explore the effect of taxonomic uncertainty on modelled pest risk. The case study involved a high profile plant pathogen Puccinia psidii, which causes a major disease of plants within the Myrtaceae family. Consequently, the results and recommendations may attract a wide range of interest in the biosecurity and pest risk modelling communities. We found the study by Elith et al. (2013) included a number of methodological issues that limit some of the specific and general conclusions reached in the paper. We discuss these issues and the ensuing implications for biosecurity management. We also draw attention to the need for pest risk modellers and biosecurity managers to find ways to communicate more effectively. We urge modellers and managers alike to develop a better understanding of the challenges and limitations of modelling species potential distributions across novel climates, and to be able to appreciate the meanings and limitations of models framed in different ways.

Biosecurity, eucalyptus rust, extrapolation, guava rust, MaxEnt, methods, myrtle rust, niche modelling, Puccinia psidii , species distribution modelling, Uredo rangelii

In a recent issue of Australasian Plant Pathology, Elith et al. (2013) presented modelling results on the high profile plant pathogen Puccinia psidii, which causes a major disease of plants within the Myrtaceae family. The disease is commonly referred to as guava or eucalyptus rust. This case study was used to explore the effects of taxonomic uncertainty on pest risk assessments for the Australian continent, and it was concluded that the estimation of pest risk can be highly sensitive to different taxonomic delimitations of the organism under investigation. The model variations either considered Uredo rangelii (commonly referred to as myrtle rust) a different species to P. psidii sensu stricto (s.s.), or a minor morphological variant within the broader species complex P. psidii sensu lato (s.l.) (Simpson et al. 2006; Carnegie and Cooper 2012). Our aims with this comment are to contribute constructively to modelling best practice for pest risk assessments, and to highlight our concerns with using the results presented in Elith et al. (2013) as a basis for estimating geographical pest risks from P. psidii. We first highlight what we consider to be a number of methodological limitations in the models of Elith et al. (2013), and give context to the attribution of causation that follows. Considering the importance of, and interest in this pathogen in Australia and worldwide, we also discuss the implications of the results and conclusions presented in Elith et al. (2013) for biosecurity management.

For at least one hundred years, several names had been applied to different populations of the rust fungus found on Myrtaceae in South America, now all regarded as P. psidii s.l. About thirty years ago, it was found that two different types of urediniospores are present in voucher specimens of P. psidii (Walker 1983). In some cases, the urediniospores are echinulate all over their surface; in the others, there is a smooth patch (called a tonsure) free of spines on the spores. These observations led to the proposition that P. psidii may be a complex composed of different morphological, and perhaps physiological, forms. In 2006, Simpson et al. (2006) named the entity within P. psidii with tonsured urediniospores Uredo rangelii, based only on two herbarium specimens that did not have teliospores. This new species determination is controversial, since the presence of a tonsure on urediniospores in some populations (with or without teliospores present) is the only distinct morphological variant found within P. psidii s.l. It is widely believed that more research on the morphology, host range and molecular characteristics of different populations is necessary before it can be decided if the variant with tonsured urediniospores and other putative variants should be considered as different varieties or distinct species. Until these questions are resolved, the rust on Myrtaceae can be referred to as P. psidii s.l., or the P. psidii complex.

The modelling system (MaxEnt; Phillips et al. 2006) that Elith et al. (2013) used to support the sensitivity analysis of pest risk maps is a commonly used correlative model. MaxEnt fits maximum entropy statistical models to relate geographical distribution points to spatial covariates. There are three issues worth elaborating regarding the development of the models presented in Elith et al. (2013), and consequently the interpretations drawn from them.

Covariate values and extrapolation

Figure 2 of Elith et al. (2013) indicates that the xeric and alpine regions of Australia are ‘suitable’ for U. rangelii, which, on the balance of evidence seems highly improbable. We suspect that both the cold and dry tolerance limits in the MaxEnt Uredo models are strongly influenced by a single misleading location record in north-western Argentina that was included in the training dataset. This record appears to be a coarsely geocoded point location representing the centroid of Tucuman, a small province that spans a large altitudinal (and temperature and rainfall) range at the base of the Andes. We think it is likely that the centroid falls in colder and drier conditions in the gridded climatology than where the specimen was observed. Similar improbable projections for the rest of the world (Elith et al. 2013, Online Appendix 3) are also not discussed, despite their potential implications for pest risk management outside Australia.

MaxEnt requires the modeller to specify how the model should fit covariate response functions beyond the range of conditions experienced by the training dataset (Fig. 8 in Webber et al. 2011; Webber et al. 2012). In the examples presented in Elith et al. (2013; Fig. 3, Online Appendix 5), it is clear that the “clamping” option was chosen. This choice means that the fitted response value for the most extreme covariate value was held constant when the model was extrapolated beyond the training data. Biologically, this default option in MaxEnt defies the Law of Tolerance (reviewed in Shelford 1963). Furthermore, no threshold was set to divide modelled suitability values into ‘suitable’ and ‘unsuitable’ regions. To explore the implications of these choices and the form of the response curves more generally, we constructed a conceptual sketch of how projections might have been influenced by the modelled relationships (Fig. 1). In doing so, we recognise that the MaxEnt ‘Explain’ tool is also particularly informative for this type of interrogation, should all the input data be available.

Figure 1.

Risk maps of modelled climate suitability and selected fitted response curves from Elith et al. (2013), for two models: Ured_27 (trained with 27 locations from records classified as Uredo rangelii) and Pucc_94 (trained with 94 location records for specimens classified as Puccini psidii). See Elith et al. (2013) for further methodological details. The cross-hatched areas are novel compared with the covariate background used to generate the respective models (based on MESS map analysis; Elith et al. 2010), indicating areas of model extrapolation. Maps and fitted response curves have been annotated with examples to indicate the likely impacts of the chosen clamping and model response curve constraints; a aridity b precipitation of the driest month c and d minimum temperature of the coldest month. Inland areas of the Ured_27 model are modelled unfeasibly as having areas of relatively high suitability in excessively dry climates. In c and d, the fitted quadratic response functions asymptote, rather than intersect the origin, leading to biologically unfeasible modelled climate suitability in excessively cold areas (Kriticos et al. 2013).

Developers of MaxEnt have called repeatedly for “extreme care” when extrapolating to novel climates (e.g. Elith et al. 2011), going as far as stating that extrapolation is “inherently risky” (Elith and Leathwick 2009; Elith et al. 2010). As such, like many correlative species distribution models, it has been shown to be poorly suited to addressing questions concerning novel climates, such as encountered in biosecurity and climate change applications (e.g., Kriticos and Randall 2001; Sutherst and Bourne 2009; Rodda et al. 2011; Webber et al. 2011). MESS maps (sensuElith et al. 2010) and the ExDet tool (Mesgaran et al. 2014) provide an informative way to visualise regions where the model is extrapolating, and therefore, where any model interpretations should be done extremely cautiously. Elith et al. (2013) state “we have evidence that the models do not need to extrapolate in our regions of interest”. However, the MESS maps presented in the paper (their Online Appendix 4) are based on the extent of the backgrounds (rather than the extent of the known presence points), only identify model extrapolation due to covariate range, (ignoring changes in the relationships between covariates), and use an ambiguous colour scheme to depict extrapolation in the geographical regions where biologically implausible areas are modelled as suitable within the “area of interest” (Australia; Fig. 1). These choices and the form of the fitted model response curves result in potentially unsuitable environments being modelled as suitable beyond the putative physiological limits of the species concerned (i.e. commission errors).

Elith et al. (2013) states “we use the Australian incursion of myrtle/guava rust as an example, not to argue which taxonomic interpretation, data set or model is correct, but to highlight the impact of taxonomic belief on modelled predictions”. However, there is no denying that the results presented will attract the attention of many biosecurity managers, both within Australia and worldwide. These results, therefore, fall within the context of a discussion about the degree of risks Australia faces from a recently discovered, high profile, non-native pathogen whose taxonomic status is uncertain.

Maps are potent communication devices, and two aspects of the way they have been presented in Elith et al. (2013) could have significant impacts on how they are interpreted by biosecurity practitioners. Firstly, the white regions in the maps in Figs 2 and 4 in Elith et al. (2013) may imply that a modelled suitability threshold of 0.05 was chosen, particularly as white was also used for the ‘no data’ class. Discussions with the authors confirmed that this is not the case, and the classification threshold was an arbitrary cartographic choice; hence, there were no explicit assumptions made about where biosecurity managers should ‘stop worrying’ about the threat of establishment of P. psidii s.l. These models, therefore, provide limited information for delimiting areas at risk of establishment or invasion (compare Fig. 1 of this paper with Fig. 4 of Elith et al. 2013).

Secondly, a reasonable interpretation of these maps appears to be that the extent of the geographical area of concern, at all levels of modelled suitability, is likely to be greater from U. rangelii than P. psidii s.l., despite U. rangelii distribution records having a narrower geographical range. Indeed, this is the interpretation intended by Elith et al. (2013: 49):

“Recognition ofP. psidiisensu lato (Puccinia_94) would lead managers to place lower priority on surveillance and containment in Western Australia, and to increase the focus of activities in Australia’s northern and eastern neighbours (e.g. New Caledonia). Recognition ofU. rangelii(e.g. Uredo_27) would lead managers to increase the priority for these activities in New Zealand, Tasmania and Western Australia (Fig. 2).”.

According to both set theory and ecological reasoning this result and conclusion are implausible. In ecological terms, based on the critical assumptions underlying correlative modelling methods and niche theory, the broader the range of environmental tolerances encompassed by an organism in its native range, the broader the range of conditions we might suppose it is at least capable of inhabiting in an introduced range. That is, if Area_small and Area_large refer to the environmental space occupied by taxa with a smaller or larger environmental envelope respectively, the corresponding environmental space projected to be suitable (potential niche breadth, Habitat_small and Habitat_large) should conform to the same inequality i.e., because Area_small ⊂ Area_large then Habitat_small ⊂ Habitat_large. The results presented in Elith et al. (2013) indicate the opposite, and are most surely a modelling artefact, due to the reduced information in the smaller datasets in relation to the background definition. Importantly though, the manner in which the MaxEnt model assigns probabilities to cells means that the map classes in the different maps in Fig. 2 and Fig. 4 in Elith et al. (2013) do not, strictly speaking, represent the same level of modelled suitability. That is, the values of modelled suitability represented in these maps cannot be compared directly between maps. This limitation is a subtlety which is not conveyed in the figure captions or body text of Elith et al. (2013), and we contend is not something that would be generally appreciated within the risk assessment community. Taken together, policies based on the prima facie results of the comparative modelling presented in Elith et al. (2013) would be in the wrong direction in terms of the perceived geographical extent of the threats.

From ecological theory, we expect that closely-related species (such as U. rangelii and P. psidii s.s., if they are eventually confirmed as being different species) may competitively exclude each other from otherwise suitable habitat (Hardin 1960; Mitchell and Power 2003), though for some species including an overlapping hybrid zone. This principle underlies the Enemy Release Hypothesis (Keane and Crawley 2002), which describes the direction of expected niche changes when organisms are freed of the effects of their natural enemies, including that of competitors. Evidence of this may be found in comparisons of habitat suitability from non-native and native ranges (Kriticos and Randall 2001). It also underpins cautions to pest risk analysts to consider interactions with other species (Davis et al. 1998) and the effects of land use in modifying the apparent climatic response (Bourdôt et al. 2013).

A corollary of this theory is that the inferred differences in climatic preferences based on native range sampling relate only to the realised niche (Hutchinson 1957; Soberón 2007); the fundamental niche of closely-related sympatric species may overlap strongly. Hence, for such species, taxonomic hair-splitting when selecting presence records for modelling is likely to produce conservatively-biased results. This method may be attractive for conservation-oriented modelling; however, the costs associated with errors regarding biosecurity risks are asymmetric. That is, the consequences of underestimating the potential range of an invasive alien species are that it could become established, creating perpetual unwanted pest impacts. The consequences of overestimating the risk are unwarranted resource allocation to biosecurity exclusion efforts. Ideally, both types of errors would be minimised when informing biosecurity policies, but overestimating the risks is generally to be preferred to underestimating them. However, the overestimates should result from ecologically robust modelling projections, rather than from modelling artefacts.

Based on their modelling, Elith et al. (2013) argue that significant resources should be devoted to resolving the taxonomic imbroglio of this rust fungus. We agree that it is critical generally for ecological modellers to attempt to clarify the taxonomy of the organism being modelled. This applies equally to geographical distribution data (Kriticos et al. 2003; Dupin et al. 2011; Baker et al. 2012; Bourdôt et al. 2013) as it does to other biological data about the taxon of concern. However, from a pragmatic biosecurity management perspective, attempting to narrow the perceived pest risks from P. psidii s.l. to those from U. rangelii faces two distinct challenges. Firstly, our ability to readily detect and distinguish samples of P. psidii s.l. and U. rangelii are impractical, and secondly, the management responses would be no different depending on which of these taxa are posing an invasion threat. We contend that in cases of taxonomic uncertainty, robust decision making in the biosecurity arena (sensuSimon 1991) is therefore best informed by assuming the broader taxonomic delimitation of P. psidii s.l., as has been done in Australia and elsewhere to date for surveillance and monitoring efforts. The expenditure of significant resources on taxonomic differentiation for pest risk assessments may result in considerable time and resource savings for some systems such as the Tephritids (Schutze et al. 2012), but for P. psidii s.l. it would seem of little practical value within a pest risk management context.

In order to better gauge the biosecurity implications of taxonomic uncertainty (avoiding the pitfalls in the MaxEnt example discussed here), we may be better off using true presence-only correlative models such as BIOCLIM/ANUCLIM (Booth et al. 2014) or CLIMEX Regional Match Climates (Kriticos 2012) that are likely to respond to changes in taxonomy in an ecologically conformal manner. Further, because all correlative bioclimatic models rely upon geographical data solely to infer climate suitability, they are inherently naïve with respect to competitive exclusion and the effects of enemy release. Their ability to extrapolate into novel climates is unreliable because they are not constrained to fit biologically meaningful covariate response functions (Austin 1987), and may over-fit to biased samples (Kriticos and Randall 2001; Rodda et al. 2011; Webber et al. 2011). Another option is to develop more mechanistic models using packages such as CLIMEX Compare Locations (Sutherst et al. 2007b) or NicheMapper (Kearney and Porter 2009), that have been designed specifically for producing ecologically plausible results when projected to novel climates (Sutherst and Maywald 1985; Magarey et al. 2007; Sutherst et al. 2007b; Kearney and Porter 2009; Sutherst and Bourne 2009; Webber et al. 2011; Sutherst 2013). Nonetheless, even armed with models well-suited to the task, it may be difficult or impossible to calibrate models of closely-related taxa sufficiently well to understand accurately the invasion risks they each pose (Sutherst et al. 2007a).

Pest taxonomy is a subject that continues to challenge biosecurity agencies (deWaard et al. 2010; Collins et al. 2012; San Jose et al. 2013). Pest risk modelling clearly has a potential role alongside molecular biology in unravelling such taxonomic puzzles, testing taxonomic hypotheses against biogeographical data (Thompson et al. 2011). However, to earn and maintain this role, the modelling should be reliable, and founded solidly in ecology. In the case of P. psidii s.l., it is our opinion that the modelled taxonomic sensitivity of pest risks presented in Elith et al. (2013) reflect significant modelling artefacts.

The specific issues we raise in this paper fall within the context of a relatively immature, rapidly evolving field of science, where methods are being adapted, developed and tested at such a rate that a consensus view of best practices has yet to emerge. As has been so carefully emphasised in the past by the authors of both this paper and Elith et al. (2013), modellers need to match the question being addressed with the most appropriate techniques available. Whilst this need applies generally, it is particularly pertinent for studies involving biological invasions or climate change where modellers are challenged with unstable range dynamics and projecting results for novel environments (Sutherst and Bourne 2009; Webber et al. 2011). The complexity of biological invasions and the related risk management questions mean that this matching process demands a significant level of modelling experience and expertise.

Some recent developments in computing technologies have been focused on making ecological modelling tools accessible to the masses (e.g., Graham et al. 2010). Even casual scans of ecological modelling discussion-lists reveal the alarming frequency with which scientists with little or no prior experience or training in ecological modelling grasp these techniques and apply them as a means of ‘rounding out their studies’. The dire consequences are apparent in the proliferation of poorly-founded and ecologically-implausible models appearing in the recent species distribution and niche modelling literature. Clearly, making it easier to generate models does not ensure that the models will be meaningful or fit-for-purpose. Modelling species’ current and potential distributions is a complex, difficult task, and there are few formal opportunities for learning appropriate modelling skills.

For end users there is a need to become more familiar with how various modelling choices affect the meaning and utility of the model results. Moreover, modellers have a responsibility to foster an effective understanding of these issues amongst biosecurity risk managers. For example, it is not possible to look at a map of a modelled species distribution, and to know instinctively what it means in terms of pest risks. There are many different ways in which the risk modelling problem can be framed, and the meaning of the model results changes accordingly.

There is clearly still much work to be done in this space, and we will need contributions from both the ecological modelling and biosecurity communities to achieve our goals of advancing best practice for pest risk modelling.

We thank John Walker for contextual information on the taxonomic history of P. psidii s.l. and Jane Elith and Mark Burgman for constructive discussions on their modelling decisions and interpretation, clarifying the background definitions and making available some of the data from their models to allow us to interrogate the output in greater detail and to construct Figure 1.