IntroductionTop

The fall armyworm is one of the pests that cause the greatest damage to the corn (Zea mays L.) crop, leading to a loss of about 34% of the grain production (Cruz et al., 1999). The insects feeds in the whorl of corn from its emergence to the formation of spikes, and if not controlled causes the death of the plant. In the world, corn is one of the most valuable cereal grain crop, it is used for both human and animal consumption. Corn is also used to make numerous non-food products. It can be cultivated in a large part of the world, and Brazil is a big world producer of the grain. Given its usefulness, it is extremely important to take into account its productive potential, resistance to diseases, types of soil besides the climatic conditions concerning to optimize the yield of grain crops.

In general, the pest control is made through the chemical pesticides, which besides polluting the environment also cause the increasing of the resistance of the insect. For this reason control techniques for the insect are researched. For example Sena Jr. et al. (2003) developed an algorithm for identifying damaged corn plants by the fall armyworm using digital colour images. Yu et al. (2003) studied detoxification capability in larval fat bodies and midguts, and adult abdomens of insecticide-resistant and -susceptible fall armyworms. They concluded that the insecticide resistance observed was due to multiple resistance mechanisms. Niu et al. (2013) studied the resistance of the fall armyworm to certain types of protein with regards to propose a hybrid resistant to the caterpillar. In Florida and Lousiana Yang et al. (2013) determined the susceptibility of S. frugiperda to a pyramided Bt corn hybrid containing Agrisure®Viptera™ 3111 corn.

Statistical methods based on spatial point processes provide useful tools for points randomly distributed on the plane or in space in many fields of knowledge, such as ecology, biology, epidemiology, seismology, archaeology, astronomy, and geography. However, most of the analyses are made assuming isotropy of the process without an adequate confirmation (Guan et al., 2006). This fact is often due to a lack of studies that consider the investigation of the isotropy in this type of data. In ecology, seismology, and biology, Guan et al. (2006), Mateu & Nicolis (2015) and Rajala et al. (2016) presented a few studies for anisotropy. The first one, proposed a formal test for isotropy based on the asymptotic joint normality of the sample second-order intensity function. Mateu & Nicolis (2015) proposed a wavelet-based test, and Rajala et al. (2016) developed a two-stage non-parametric method for quantifying geometric anisotropy for patterns that are compressed or stretched.

The use of spatial point pattern analysis in agriculture was considered by Spósito et al. (2007), who studied the disease of the black spot in an orange crop in Brazil. They evaluated the point pattern of the infected trees by three different methods, the binomial dispersion index, the Ripley K function and Monte Carlo test, and concluded that the disease shows a clustered behavior in the region under study. Ribeiro Jr. et al. (2009) used point processes to detect and model specific patterns in the occurrence of the thrips disease in onion. They used the Mantel test to detect the point pattern and later described the pattern by a mixed spatial Poisson process with a geostatistical random component. They also obtained maps of the prediction of the levels of susceptibility to infestation in the study area.

In this work we study the spreading of the fall armyworm along an agricultural experimental area at Paraná State, Brazil, comparing methods of anisotropy by means of the spatial point patterns. The main purpose is to examine its spatial distribution along the study area, investigating possible dominant direction of infestation of the fall armyworm, and if the pattern presents characteristics of clustering, regular or random spreading. This knowledge provides the farmer a more adequate management of his crop by considering the characteristics of the dispersion of the insect along the area in consortium with the techniques of precision agriculture. The lack of literature that incorporates spatial point processes and agriculture motivated this work.

Material and MethodsTop

The data set

The experiment was carried out in an agricultural commercial area of 27 ha located in the city of Cascavel, Paraná state, Brazil (approx. 24.95o S, 53.57o W and 650 m asl). Local soil is classified as Rhodic Hapludor (Soil Survey Staff, 2014). The climate of the region is classified as mesothermic and super humid temperate with average annual temperature around 21oC. A representation of the location of the experimental area is presented in Fig. 1a.

Figure 1. Experimental area in Cascavel, Paraná State, South of Brazil (a), and satellite image of experimental area with the polygon of the study area (b).

For the data sampling, we considered as a possible event each plant of corn attacked with the fall armyworm in the experimental area during the crop year 2015/2016. For every observed infected plant we registered the geographical coordinates as an event of interest, which were obtained using a GEOEXPLORE 3 GPS positioning system receiver in a Universal Transverse Mercator (UTM) coordinate system.

The study area (Fig. 1.) included 9360 corn plants of which 1303 were infected by the S. frugiperda. This study area was cultivated with the hybrid Pioneer 30F53YH variety, developed by the company DuPont Pioneer, under direct seeding system, with a population of 82,222 plants/ha. Along the season we applied 1095 kg of fertilizers (8% N, 20% P and 18% K). The planting of the seeds began on August 28, 2015 and ended on August 31, 2015 and was performed with line spacing of 0.45 m and between plants of approximately 0.33 m. The sampling was made on October 15, 2015, when the corn crop was in vegetative growth stage.

Statistical analysis of anisotropy

A spatial point pattern is a stochastic mechanism that generates a countable set of events xi in the plane (Diggle, 2003). Events can be scattered in the region of the plane in three forms of relationship, namely completely spatial random (CSR) as represented in Fig. 2a, regular in Fig. 2b, and clustered in Fig. 2c. In spatial point analysis, the main objective is to identify possible properties and forms of relationship of the phenomenon under study taking into account its spatial sites, through geographical coordinates.

Figure 2. Simulations of Poisson spatial point patterns, random (a), regular (b) and cluster (c). Estimated intensity function for random (d), regular (e) and cluster simulated pattern (f).

A spatial point process is considered to be stationary under translations if the distribution is unchanged when the origin of the index set is translated (Ripley, 1981; Diggle, 2003; Illian et al., 2008), and is said to be isotropic if its distribution is invariant under rotations about the origin. Otherwise the process can be said to be anisotropic. The assumption of stationarity simplifies drastically the statistics of point patterns (Baddeley et al., 2006) because the chance of observing some point configuration at a specific location is independent of the location. Usually, the assumption of isotropy is made by simpler interpretation and ease of analysis (Guan et al., 2006). In other words, stationarity implies that it is possible to estimate some properties of the process from a single realization on the region A, because these properties are the same in different, but geometrically similar, sub regions of A. Isotropy means that they are invariant with respect to directions (Mateu & Nicolis, 2015). In the following, we assume that the process N is stationary, but not necessarily isotropic. The simplest example of a stationary and isotropic point process is the homogeneous Poisson point process.

Figs. 2(a), 2(b) and 2(c) show simulated CSR, regular and cluster spatial point patterns, respectively. The difference of them is the distance between the points. For the regular pattern the distances between neighbor points appear to be constant, while for the cluster pattern there are groups of points over the region. In the CSR pattern the distance between the points is random.

Consider N, a two-dimensional spatial point process, where each event is given by x ∈ R² . Let A be the region where the points are taken, denote as |A| the area of this region, and as N(A) the random number of points in A. The intensity function may be described by the first order property, given as

where E[.] is the expected value; N(.) is the number of events in the plan region; dx is a small region containing the event x; |dx| is the area of the region dx. For a stationary process λ(x) = λ , where λ is the average number of events per unit area (Daley & Vere-Jones, 2002).

The estimation of the first-order property, or intensity of the point pattern, is usually performed through kernel functions (Baddeley et al., 2006, 2015),

where k(.) is the kernel function with ∫ k(u)du = 1, and h>0 is the smoothing parameter. The intensity estimate at the point x depends only on the spatial relationship between x and the elements of the sample, xi (i = 1, …, N(A)), quantified by the metric embedded in the kernel function. In order to estimate λ(x) by using the kernel function as described in Eq. [2], we used the kernel2d function of the package splancs (Rowlingson & Diggle, 2016) represented by Figs. 2d-2f. The kernel intensity function estimation using quartic smoothing parameter, h=0.2, was performed for each simulated pattern and it is shown in Figs. 2d, 2e and 2f corresponding to the estimation of the first-order property for each type of simulated point pattern. It is worth to note that the intensity seems almost constant for the simulated regular pattern as in Fig. 2, and for the simulated cluster pattern there are regions with different values of intensity as in Fig. 2f.

The second-order properties of a point process can be observed using the covariance between events of two areas. The second-order intensity provides a way of describing the relationship between pairs of points and is defined as

where dx is a small region containing the event x; dy is a small region containing the event y and |dx| and |dy| are the areas of regions dx and dy, respectively. For a stationary process, λ2 (x,y) = λ2 (x-y), and for a stationary, isotropic process, λ2 (x-y) = λ2 (r) where r =||x-y|| is the Euclidean distance (Daley & Vere-Jones, 2002).

Ripley (1976, 1977) presented an alternative characterization of the second-order properties (Eq. [3]) of a stationary and isotropic process, given by

The term λK(r) is the mean number of points other than the typical point in a ball of radius r centered at the typical point. The K :[0, ∞ ) → [0, ∞) function is estimated as

where dij is the distance between the i-th and j-th points, I(x) is the indicator function with the value 1 if x is true and 0 otherwise, and λ is the estimate of Eq. [1].

This function has shown its importance as a spatial dependency summary since it allows quantifying spatial dependence between different regions of the process, and is also known as reduced second-order analysis. It has been widely used in two-dimensional point processes because it enables the detection of the point pattern in different distance scales simultaneously and the observed point pattern can be compared to known stochastic process models for different point configurations (Diggle, 2003).

The shape of K(r) relative to that of the CSR provides valuable information on the point process distribution (Illian et al., 2008); for a CSR process we have K(r) = pr2. For a cluster process the points have more neighbors than expected under CSR, and therefore the estimates for K(r) will be greater than pr2. On the other hand, for a regular pattern, the estimates of K(r) will be less than pr2 (Diggle, 2003). Fig. 3 represents the K function (Eq. [4]) and the empirical envelope of 95% for each pattern simulated in Figs. 2a, 2b and 2c, respectively.

Figure 3. Ripley K function for the simulated point patterns: random (a), regular (b) and cluster (c). The red dotted line represents the theoretical CSR case K(r)=Πr2and the grey region is the empirical envelope of 95%. The continuous black line represents the observed data.

A directional counterpart of Ripley K function was proposed by Rajala et al. (2016) as a non-parametric methodology to study anisotropy in spatial point patterns. Instead of counting all data points falling inside a circle of radius r centered at a data point, we could replace the circle by another geometrical shape (Baddeley et al., 2015). For a unit vector u, the conical directional K-function is defined as

where λ is the intensity, and Cu(r,Θ) denotes a double cone in the direction u with an slant height of length r and an apex angle of size 2θ centered in 0 (Redenbach et al., 2009).

The function is obtained by counting how many pairs of events in the pattern have both their vector of difference angle less than θ radians from direction u and difference vector length less than range r. The procedure consists in comparing each of the directions that is defined by the researcher with the function under CSR. In order to do this, it is necessary to assume that anisotropy is geometrical, i.e. a result of compression or stretching and modeled by a linear transformation comprised of scaling and rotation (Rajala et al., 2016).

Mateu & Nicolis (2015) developed an anisotropy wavelet-based test by using the empirical logarithm of the directional scalogram, assuming that under the null hypothesis of isotropy a random point pattern is expected having the same value of the directional sca-logram for all possible directions. Let f(x), x ∈ R2, be a function. The continuous directional wavelet transform for scale a and orientation θ is given by

where the over line denotes complex conjugate. There are many different directional wavelets ψ a,b (x, θ) proposed in literature, see for example Neupauer & Powell (2005).

Kumar (1995) proposed the function

with regard to identify the behavior of the process in different directions. η (a, θ)characterizes the distribution of the energy at different scales and directions. The directional scalogram is obtained from the variance of wavelet coefficients for each scale and each direction, given by

Mateu & Nicolis (2015) applied the fully anisotropic Morlet wavelet transforms given in Neupauer & Powell (2005) to the intensity function λ (x), estimated by the box counting method, where x denotes locations. Let m be the possible directions θ _t ∈ (0,π] with i=1, …, m; L possible scales aj with j=1, …, L, and N(A) spatial points bk, k=1, …, N(A). The directional wavelet transforms

are conducted for a range of scales and orientations at all positions in the domain of λ (x). Denote by S_{θi, aj} the variance of the corresponding wavelet coefficients for a particular direction θi and scale aj,

Monte Carlo simulation determines the significance of the directional wavelet coefficients. By assuming that the null hypothesis is isotropy, an isotropic random point pattern should be expected to have the same value of the directional scalogram for all possible directions. Nicolis et al. (2010) considered the logarithm of the variance of the wavelet coefficients over all scales as test statistic,

An alternative method to study the isotropic behavior of the data set is the histogram of the k-th nearest neighbor angle. For each event in the point pattern, we calculate the distance to the k-th nearest neighbor and the azimuth angle a between this pair of events where a α ∈ [-π, π]. The isotropic behavior of the data set is done by comparison with the theoretical value 1/2p where peaks above this theoretical value indicate a possible dominant angle direction of anisotropy. The histograms can be used to investigate in which nearest neighbor the process reaches isotropy, that is, when there are no more peaks in the histogram. We used the knnangle function of the Kdirectional package (Rajala, 2017) to build the histograms.

To confirm the form of relationship of the point pattern we calculated the quadrat counting test (Illian et al., 2008). It is a test of CSR that uses the χ 2 statistic based on quadrat counts where the study area A is divided into sub regions called quadrats (A1, A2, …, Aq) of equal area. The test counts the number of points that fall in each quadrat nj = N(A n Aj) for j=1, ..., q. Under the null hypothesis of CSR, the nj are independent and identically distributed Poisson random variables. The test statistic is given by

where s2 is the sample variance. If the distribution is Poisson, the variance is exactly the mean under CSR s2 / = 1 . Moreover, if s2 / < 1 then there is too little variation among quadrat counts, suggesting regularity rather than randomness. Similarly, if s2 / >1 there is too much variation among counts, suggesting clustering.

For the data analysis we used the software R version 3.3.3 (R Team, 2017) and the R packages spatstat version 1.47-0 (Baddeley & Turner, 2005) and Kdirectional version 1.01 (Rajala, 2017).