Conservation Strategy for Palm Groves: Optimal Chemical Control Model for Red Palm Weevil, Rhynchophorus ferrugineus

Abstract

Rhynchophorus ferrugineus (Olivier, 1790) is an invasive pest species that constitutes one of the most important problems around the Mediterranean region and has been responsible for the loss of over 100,000 palm trees with an estimated annual cost of EUR several hundred million since its introduction into Europe. Methodological approaches of conservation ecology, such as multidisciplinary modelling, also apply in the management of cultural landscapes concerning ornamental plants, such as palm trees of the area. In this paper, we propose a dynamic model for the control of the red palm weevil, contributing in this way to the sustainability of an existing cultural landscape. The primary data set collected is a sample from the density-time function of a two-cohort pest population. This data set suggests a bimodal analytic description. If, from this data set, we calculate a sample from the accumulated density-time function (the integral of the density-time function), it displays a double sigmoid function (with two inflections). A good candidate for the analytical description of the latter is the sum of two logistic functions. As for the dynamic description of the process, a two-dimensional system of differential equations can be obtained, where the solution’s second component provides the analytical description of the original density-time function for the two-cohort population. Since the two-cohort waves appear in all three cycle stages, this reasoning applies to the subpopulations of larvae, pupae and adults. The model fitting is always performed using the SimFit package. On this basis, a mathematical model is proposed, which is sufficiently versatile to be of help in the control of this pest species in other geographical areas.

1. Introduction

The red palm weevil, Rhynchophorus ferrugineus, hereinafter referred to as RPW (Olivier, 1790) (Col.: Dryophthoridae), is a species native to Asia and Polynesia, and is characterized as an invasive pest of great severity and economic importance worldwide [1,2]. This species has expanded from its area of origin to different continents. It has been recorded in 28 countries in Asia, six in Africa, one in North America, two in Central America and the Caribbean, 14 in Europe and five in Oceania, causing serious problems in coconut trees, both in crop cultivation and ornamental use [1,3].

In Spain, RPW has spread widely over the Mediterranean coast, the interior of the peninsula and the islands, causing significant ecological and economic damage in areas such as the “Palmeral de Elche” (in the Valencian Community), declared a World Heritage Site by UNESCO, and in the Phoenix canariensis forest on the Canary Islands, which is a natural heritage site and represents the most relevant source of genetic variation for this palm tree [1,4,5].

Many aspects of RPW biology and ecology have been widely studied and are reported in several publications (e.g., [1,6,7,8,9]). However, the cold survival strategy adopted by this species, native to warm geographic areas, remains poorly known [10].

RPW infestations threaten agriculture (e.g., date plantations), affect the aesthetics of infested palm landscapes and uncontrolled urban populations (e.g., UNESCO World Heritage palm groves) and threaten natural areas (e.g., palm oases and endemic palm) [1,8,9,11,12].

In Spain during the 2005–2009 period, between 47,000 and 50,000 infested palm trees were destroyed; this involved containment/eradication costs of EUR 44.5 million, particularly in the Valencian Community where professional and family managed nurseries were affected [11,13,14]. In this regard, the FCEC [11] considered that infested (but not destroyed) palm trees might result in losses of between EUR 96 and 288 million.

Insecticide applications are, at present, the most effective method for protecting palms from attack by palm weevils [8]. Larval control, mainly by injection of insecticide into the trunk, carried out during the early stages of infestation, sometimes accompanied by application of spray for adult control, is currently the most employed (e.g., [3,8,15,16,17,18,19]). This is usually accompanied by monitoring and mass trapping of adults using semiochemicals (e.g., [3,9,20]). It should be noted that other control methods have been studied and developed, e.g., biological control (with entomophages and entomopathogens) and genetic methods (sterile backcrosses, sterile insect technique and chemosterilization), including two or more of the above control measures in integrated pest management programmes (e.g., [3,8,20]). The severity of R. ferrugineus and the problems associated with its control, including the difficulty of detecting infested palm trees early, the risk of mass captures in non-infested areas and the loss of control strategy effectiveness [9], motivate the search for new tools that increase the bioecological knowledge for dealing with RPW populations. In this regard, mathematical models can offer a holistic approach to understanding the biological behaviour of species, and they are an instrument for predicting and/or simulating complex systems that include the population dynamics of the pest and its control strategies [21,22,23].

In the plant protection field, mathematical models have focused on using sigmoid functions, such as Verhulst–Pearl or logistics, Gompertz, Weibull, Richards and monomolecular, as well as those models that combine them [24,25,26,27]. These functions have also been integrated to optimize the use of the resources and time spent on combating pests, through the formulation of optimal control problems aimed at keeping a pest population below the economic damage level as well as reducing the need to implement control strategies [21,22,23,28,29].

The existing literature on mathematical models related to the ecological aspects of R. ferrugineus is scarce, with only the work of Rossini et al. [30], who proposed a linear model for RPW development as a function of temperature. However, Hansen et al. [31] adapted a logistic model for another Curculionidae species, Dendroctonus rufipennis (Kirby), based on temperature as a diapause prediction factor of the insect’s fourth larval stage. It should be noted that the larval diapause in Curculionidae has also been recorded for Curculio elephas Gyllenhal, C. sikkimensis (Heller) and Smicronyx fulvus LeConte, where the larvae were subjected to temperatures below 15 °C [32,33,34,35].

In the R. ferrugineus larvae, Martín and Cabello [1] recorded a slowdown in the development of their last stages when the temperature dropped and remained in the 10–15 °C range for a period of more than 80 days, arguing that this phenomenon was due to a larval diapause, because the delay in the larval development and metabolic activity was related to changes in temperature [36,37].

The importance of RPW within palm production systems and/or ecosystems, and the need to deepen our understanding of the population behaviour of its stages, is fundamental to the development of mathematical models that can describe the population dynamics of the larvae, pupae and adults of R. ferrugineus, and to propose an optimal solution to control their larvae. This was the objective of this work.

The resulting modelling–methodological development also contributes to the sustainable management of the involved cultural landscape of palm groves and ornamental palm trees.

The objectives of this work were to answer two initial issues:

(1)

Is it possible to fix a mathematical model for the biology of RPW in conditions of cold winter climates?

(2)

Once the above model has been fixed, is it possible to use an optimal control model to simulate the optimum time for application of chemical larval control?

2. Materials and Methods

2.1. Biological Data Used

The trial was carried out in the “Cortijo del Olivar”, Municipal District of Almuñecar, Granada, Spain, 50 m above sea level (13.748219 N, 3.690619 W). It was in a room, communicated with the outside, located inside a country building; it was carried out between July and September of the following year. The trial started with a total of 70 RPW eggs, less than 24 h old, which were obtained from the laboratory rearing populations (Entomology Laboratory, University of Almeria, Almeria, Spain). The food substrate for larvae and adults was sugar cane (cultivar: “Algarrobeña”; traditionally grown in southern Spain, Almuñecar). Specimen handling was carried out following the methodology described in [38,39]. After development, 13 pairs were used for the evaluation of adult parameters. Data recorded were number of days for egg, larval and pupal stages, as well as mortality at these stages; adult longevity and fecundity were also recorded. Temperature was recorded daily by thermograph.

2.2. Mathematical Model

The primary data set collected is a sample from the density-time function of a two-cohort pest population. This data set suggests a bimodal analytic description.

Originally, the single logistic function is:

$$\frac{K}{1+\left(\frac{K}{x_0}-1\right)e^{-rt}}$$

(1)

where x₀ is the initial population; r: the logistic growth rate or steepness of the curve; K: the curve’s saturation value; in function (4) below, i subscripts (=1, 2) indicate the first or second cohort of the population.

Function (1) is the solution of the logistic differential equation:

$$x^{\prime }=rx\left(1-\frac{x}{K}\right)$$

with the initial condition $x\left(0\right)=x_0$. Function (1) is sigmoid; its derivative, with appropriate constants $a_0, a_1, a_2>0$, is

$$y_1\left(t\right)=\frac{a_0{e}^{-a_{1}t}}{{\left(1+a_2{e}^{-a_{1}t}\right)}^2},$$

(2)

displaying a unimodal curve. The sample from the density-time function of a two-cohort pest population suggests a bimodal analytic description, and for this, the sum of two functions of type (2) is a promising candidate [25]:

$$y_2\left(t\right)=\frac{a_0{e}^{-a_{1}t}}{{\left(1+a_2{e}^{-a_{1}t}\right)}^{2 }}+\frac{a_3{e}^{-a_{4}t}}{{\left(1+a_5{e}^{-a_{4}t}\right)}^2},$$

(3)

where $a_0, a_1, a_2, a_3, a_4,a_5>0$ are the fitting parameters.

If we calculate a sample from the accumulated density-time function (the integral of the density-time function) from the data set of the density-time function, it shows a double sigmoid function (with two inflections). From the above, we conclude that a good candidate for the analytic description of the latter is the sum of two logistic functions (1):

$$\frac{K_1}{1+\left(\frac{K_1}{x_{01}}-1\right)e^{-r_{1}t}}+\frac{K_2}{1+\left(\frac{K_2}{x_{02}}-1\right)e^{-r_{2}t}}$$

(4)

Here, $K_i, r_i,x_{0i}>0$ (I = 1,2) are the fitting parameters.

It is worth emphasizing that, if we consider a two-cohort pest population of Rhynchophorus, the above reasoning can be applied to all the developmental stages: the larvae, pupae and adults.

2.3. Dynamic Model of Optimal Chemical Control

Considering any given stage, from functions (2) and (3), we easily obtain the following system of differential equations, where $y_1$ is the density-time function of a single cohort, and $y_2$ is the density-time function of the two-cohort population:

$$y_1^{\prime }=a_1{y}_1\left(1-\frac{2}{1+a_2{e}^{-a_{1}t}}\right)$$

(5)

$$y_2^{\prime }=a_1{y}_1\left(1-\frac{2}{1+a_2{e}^{-a_{1}t}}\right)+a_4\left(y_2-y_1\right)\left(1-\frac{2}{1+a_5{e}^{-a_{4}t}}\right)$$

(6)

Although the dynamic model (5)–(6) can be adapted to any stage, for the chemical control of Rhynchophorus, the system (5)–(6) is used for a two-cohort larval population [1]. An operational time interval is fixed [0,T], and the set of admissible controls ${\mathrm{U}}_{\mathsf{\epsilon }}\left[0,\mathrm{T}\right]$ consisting of piecewise constant functions over a fixed uniform division of interval [0,T] with $0\le u\left(\mathrm{t}\right)\le \mathsf{\epsilon }$ is defined. Constants ${\mu }_1$,${\mu }_2>0$ express the efficiency of the chemical on the respective populations. Then, from the (5)–(6) dynamics, we obtain the following control system:

$$y_1^{\prime }=a_1{y}_1\left(1-\frac{2}{1+a_2{e}^{-a_{1}t}}\right)-{\mu }_{1}u y_1$$

(7)

$$y_2^{\prime }=a_1{y}_1\left(1-\frac{2}{1+a_2{e}^{-a_{1}t}}\right)-{\mu }_{1}u{y}_1+a_4\left(y_2-y_1\right)\left(1-\frac{2}{1+a_5{e}^{-a_{4}t}}\right)-{\mu }_{2}u{y}_2$$

(8)

Now, on the basis of this control system, we set up the following optimal control problem: suppose that the damage caused by larval density $y_2$ in unit time is $\alpha y_2$; then, the total damage during the time interval [0,T] is ${\int }_0^T\alpha y_2\left(t\right)dt$. If this damage is remains below a boundary, H, for the palm tree to survive, and we wish to minimize the cost of the chemical control proportional to ${\int }_0^{T}u\left(t\right)dt$, we obtain the following constrained optimal control problem:

$$\left.\begin{array}{c}\Psi \left(u\right)={{\displaystyle \int }}_0^{T}u\left(t\right)dt\to min,\\ u\in U_{\epsilon }\left[0,T\right]\\ \mathrm{Equations} \left(5\right)-\left(6\right) \mathrm{hold}\\ {{\displaystyle \int }}_0^T\alpha y_2\left(t\right)dt\le H\\ y_1\left(0\right), y_1\left(0 \right) \mathrm{are} \mathrm{given.}\end{array}\right\}$$

(9)

2.4. Software Used

The above-mentioned mathematical models were fitted to the data using the TableCurve 2d v 5.0 [40] and SimFit v 6.9.8 [41] software packages.

3. Results

The results found were used to carry out the adjustment to the proposed mathematical model that explains the population dynamic of the insect’s life-cycle stages (Section 3.1). Subsequently, with this model, a control variable was used that would allow us to set the time for the chemical control of the pest in order to obtain an optimal result.

3.1. Results of Model Fitting

The results show that a portion of the larval population of R. ferrugineus developed normally (relative maximum in tp₁), while another portion of the population lengthened its development up to 370 days (

Table 1. Values of fitting parameters of model (3) for the all stages of pest population of Rhynchophorus.

Stage	Fitting Parameters (*)						Statistical Parmeters
	a0	a1	a2	a3	a4	a5	d.f.	R2 adj.	p
Larvae	673.739	0.04410	3.56691	2117.26	0.01967	80.0043	38	0.97654	<0.01
Pupa	74980.4	0.05752	792.265	10500.9	0.02003	115.495	40	0.97744	<0.01
Adult	28372.7	0.04231	485.713	218719	0.02546	2901.76	40	0.95999	<0.01

$$(*)a_0=\left(\frac{K_1}{X_{01}}-1\right)\cdot K_1\cdot r_1;a_1=r_1;a_2=\left(\frac{K_1}{X_{01}}-1\right)$$

$$a_3=\left(\frac{K_2}{X_{02}}-1\right)\cdot K_2\cdot r_2;a_4=r_2;a_5=\left(\frac{K_2}{X_{02}}-1\right)$$

and Figure 1). In the pupal and adult stages, two relative maxima (tp₁ and tp₃) are shown with a distance between them of approximately 160 days.

We note that the R² value indicates that the goodness of fit is better (higher) for the double logistic function fitted to the accumulated data (0.99951) than for the fitting dynamics (5)–(6) (0.90319). Similar results were obtained for the rest of the analysed data. This observation justifies the use of a sigmoid function fitted to the accumulated data, smoother than the original density data, which displayed more random oscillations.

Subsequently, double logistic curves were also fitted to the accumulated density data (see

Table 2. Values of fitting parameters of double logistic functions for pest populations of Rhynchophorus.

Stage	Term (i)	Fitting Parameters			Statistical Parameters
		Ki	x0i	ri	d.f.	R2 adj.	p
Larvae	1	2951.43	214.53	0.06771	37	0.99764	<0.01
	2	1673.39	57.6115	0.01661
Pupa	1	1810.51	4.8024	0.05484	40	0.99872	<0.01
	2	4621.58	47.1067	0.01889
Adult	1	1342.24	1.2492	0.04789	43	0.99951	<0.01
	2	3031.58	1.32482	0.02461

and Figure 2). In the latter, we also indicate the time moments t_k, where the minima and maxima of the densities are attained, according to the inflexions of the accumulated densities.

3.2. Solution of the Optimal Control Problem

For a palm tree to survive, our objective is to keep the effect of the larval population under a boundary, H, at minimum cost. Using the Table 1 fitting parameters, we solved the optimal control problem (9) with H = 2000 units of biomass, T = 450, ${\mu }_1=0.95$, ${\mu }_2=0.5, \alpha =0.7, \epsilon =0.1, y_1\left(0\right)=13, y_2(0$) = 14.

We solved the optimal control problem applying the toolbox developed for MatLab by Banga et al. [42], the results of which are shown in Figure 3.

4. Discussion

Aplasticity is well understood when it takes the form of responding to predictive environmental signals received during the individual’s lifespan. However, the so-called mixture of intragenotypic strategies has been neglected, perhaps due to the general assumption that it must be non-adaptive [43]. In this regard, several species of insects, especially from the Coleoptera order, may present two types of larval diapause (simple or prolonged), leading to diverse cohorts that differ in their life cycles, especially in terms of duration. The causes of what Ushatinskyaya [44] called “superdiapause” are poorly studied and not understood in detail. In our case, it seems clear from the results that RPW belongs to this group of insects.

The results for the first species, R. ferrugineus, showed that a portion of the larval population develop normally, but another portion could lengthen their development up to 370 days (Figure 1); this is demonstrated by the two relative maxima found in the pupal and adult stages. The above statement might indicate that the diapause exists as a survival mechanism. This phenomenon was later verified in other trials [1]. The presence or lack of a diapause in this species is crucial for the early detection and control of the pest, as well as for border inspections, since detection is carried out by ultrasound—if the larvae are in diapause, they are not feeding, and, therefore, do not make a noise with their jaws that can be detected [45]. Similarly, the results can be used to develop more complex predictive models to track the dynamics and expansion of this pest species.

The results of the logistic model developed for the population dynamics of the RPW larvae show a delay of larval development that would relate to the diapause manifesting as a survival mechanism. Martín and Cabello [1] noted that the diapause of RPW larvae was stimulated by a decline in temperature, of values between 10 °C and 15 °C, during the larval growth stage. Additionally, they indicated that this phenomenon promoted the division of the larval population into two groups: the first formed before the temperature drop occurred, in which the larvae continued their normal development (such as that recorded in tp₁; Figure 1a), while the second group delayed its development by 170 days (similar to that shown from tp₂ onwards; Figure 1a), with the consequent delay in pupal and adult morphogenesis (as in the relative maxima, tp₁ and tp₃, in Figure 1b,c, respectively).

The larval development behaviour of R. ferrugineus may represents a survival strategy that offers the maximum development opportunity in a univoltine cycle, as one group of larvae is able to accelerate their development to promote adult formation during the spring, while the other can slow down its development, influenced by the low winter temperatures, allowing adults to form that emerge during the autumn. This behaviour was also recorded by Hansen et al. [31] for the curculionid D. rufipennis, whose larval diapause was mathematically modelled in a logistic function that considers the thermal threshold as a predicting factor [31]. Other Curculionidae may also have a larval diapause in response to the drop in temperature, as is the case with Curculio elephas Gyllenhal [32,33], C. sikkimensis (Heller) [34] and Smicronyx fulvus LeConte [35].

Whether or not the RPW larvae are in diapause is crucial for early detection and pest control because most RPW inspection systems include equipment capable of detecting ultrasound produced by the larvae during feeding [46,47]. Therefore, detecting no noise will not mean the absence of larvae inside the palm stems, since diapause larvae do not feed [45].

The mathematical models developed for the larval, pupal and adult stages of RPW constitute a valuable tool for bioecological and/or prediction studies at the population level of this pest in various Spanish zones where the ambient temperature varies over the year, and can even stay below 15 °C during the winter months. At this time, the larvae can slow their development as a survival mechanism, as recorded by Martín and Cabello [1] in the municipality of Almuñecar (Granada, Spain) from the month of November to March, in which the temperature, both of the air and inside the palm tree, remained between 10 °C and 15 °C, with fluctuations within the palm tree of around 10 °C for more than 80 days.

In this regard, mathematical models can be used in various parts of Spain where palm trees are grown, such as the Andalusian coastal zones of the Mediterranean and the coastline from the Ebro Delta to Gibraltar, the Valle bajo of the Guadalquivir and Guadiana, as well as in coastal areas of the Rias Bajas in Galicia, where the ambient temperature can range from 10 °C to 12.5 °C during the month of January [48], also having an impact on the internal temperature of the palm trees.

Moreover, the optimal control model for RPW larvae reveals that the application of a chemical control measure can significantly reduce the initial larval population peak over the first 50 days. The action of chemical control also produces a significant decline in the second larval maximum, which begins to form from about day 100 in the uncontrolled population.

In turn, the rise in the chemical control trajectory seen in Figure 3, which begins after day 100, corresponds to the need to increase the use of insecticide to keep the palm biomass threshold below the economic damage level, since the larval development continues in this period, having been delayed by the temperature decline experienced by the population without chemical control. Similarly, this increase in the control trajectory can be related to the formation of RPW pupae and/or adults, which require greater control.

The logistic models developed for the population dynamics of the R. ferrugineus stages, as well as the optimal control model proposed, are instruments that can be incorporated into integrated pest control programmes, with the aim of improving decision making and reducing the cost of managing the RPW population. Similarly, these tools can be included within more complex prediction models, allowing us to understand RPW dynamics and expansion.

In this regard, it should be noted that this mathematical model has several advantages:

(A)

The model is sufficiently simple and flexible to allow its extension to this pest species in other geographical areas.

(B)

It saves some of the costs of experimentation on the effectiveness of the chemical control of RPW, as simulations can be carried out before expensive field trials, providing valuable information in advance.

(C)

With the introduction of optimisation into the model, the model could be used to calculate the ideal timing of chemical pest control, depending on climatic conditions and/or stage of infestation, with potential savings in treatment costs.

In summary, the results from the logistic models, for the pest species as well as for the different aspects of their biology/ecology, suggest that they could be an excellent tool to use in integrated pest management (IPM) programmes. This may also corroborate the results found regarding the dynamics of crop diseases [24,25], as their use can be generalized to crop protection.