INTRODUCTION
Herbicide resistance is the main problem in weed control around the globe and one of the greatest threats of modern agriculture. Accordingly, studies on herbicide resistance are at the forefront of current research in weed science (Neve et al., 2009; Panozzo et al., 2015a).
An early identification of the evolution of herbicide resistance in the field can help, if acted upon, to prevent, manage, mitigate and delay the spread of herbicide resistance (Burgos et al., 2013; Babineau et al., 2017), and in designing and implementing strategies to handle weed infestation (Panozzo et al., 2015a). An early identification of the evolution of herbicide resistance in the field can help, if acted upon, to prevent and delay the spread of resistance to neighboring weed populations.
The spatial scale on which herbicide resistant genes evolve is still poorly understood; in fact, the majority of studies on resistance only consider one or a few weed populations (Menchari et al., 2006; Schulz et al., 2014; Matzenbacher et al., 2015; Babineau et al., 2017). To handle new weed outbreaks that are resistant to herbicides, or invasive and harmful weeds, farmers need spatial information that describes the location and the propagation, within and across fields, of the problematic weeds (Cardina et al., 1997). Frequently updated and rapidly circulating maps of resistant biotypes are important in managing and tracing the evolution of herbicide resistance. Large scale maps that show the dissemination of the resistance at a continental and world level are already available in herbicide resistance global databases (Heap, 2019), nonetheless, more detailed maps are required at a national, regional and local level (Panozzo et al., 2015a).
The Echinochloa species are prone to developing resistance (Valverde et al., 2000). Echinochloa colona (L.) Link, is a hexaploid or tetraploid species (Yabuno, 1962), predominantly autogamous, which reproduces primarily through seeds but is able to reproduce through nodes (Masood et al., 2016). It is considered one of the most difficult grass weeds to control in summer and in vegetable crops in more than 60 countries (Masood et al., 2016). It is a dominant weed that is almost ubiquitous in rice crops (Valverde et al., 2000). High densities of E. colona have been reported in 24 countries with direct seeding rice crops, 12 with irrigated rice crops and in transplanted rice (Masood et al., 2016).
In rice crops, simple and multiple resistance to E. colona has been reported for the majority of modes of action (MOA), including photosystem II inhibitors (propanil), acetyl-CoA carboxylase (ACC) inhibitors (fenoxaprop, fluazifop, haloxyfop, cyhalofop-butyl), acetolactate synthase (ALS) inhibitors (penoxsulam, bispyribac-sodium, azimsulfuron, imazapyr, imazethapyr), 5-enolpyruvylshikimate-3-phosphate (EPSP) synthase inhibitors (glyphosate) and auxinic herbicides (quinclorac) (Heap, 2019). In Colombia, there are reports of resistance to propanil (Fischer et al., 1993), quinclorac (Valverde and Itoh, 2001) and penoxsulam in rice fields in Tolima (Carranza and Plaza, 2015).
Colombian farmers frequently encounter difficulties in controlling E. colona. Additionally, there is an overdependence of herbicides to control weeds and intensive planting of rice. These facts and the absence of scientific work focusing on resistance in Colombia suggest that there are unreported cases of herbicide resistance of E. colona.
The aforementioned aspects support the objectives of the current investigation: i) to determine the resistance of E. colona to the main post-emergent herbicides commonly used in rice crops, and ii) to map the spatial distribution of the resistance in the Saldaña and Purificación municipalities in Tolima.
MATERIALS AND METHODS
Plant material
In the second semester of 2016, we collected seeds from mature Echinochloa colona (L.) Link panicles in rice (Oryza sativa L.) fields nearing harvest. Each population was defined as the bulk of seeds collected in a single field. We collected 26 populations (labeled 42-67 - the labels number a larger set of populations). Each population was georeferenced and management condition, field history, weed management, and herbicides were recorded, through a structured interview. Additionally, two reference populations, previously characterized by Carranza and Plaza (2015), were also included: TA5 resistant to acetolactate synthase (ALS)- and acetyl-CoA carboxylase (ACC)-inhibitors herbicides and quinclorac, and VM susceptible to such herbicides. The seed of three of the collected populations (44, 54 and 55) did not germinate, so they were excluded from the study.
Sampling the populations of Echinochloa colona
We carried out the sampling in Saldaña (3°56’5” N, 75°1’13” W; 300 m a.s.l.) and Purificación (3°51’31” N, 74°55’52” W; 310 m a.s.l.) municipalities in Tolima state (approximate area 66.56 km2). We used a grid of 26 squares, each of 2.56 km2, to guarantee at least one square for every thousand cultivated hectares (22 000 ha at the time of sampling). The grid was designed using the geographic information system (GIS) ArcGIS 10.4.1 trial program (Environmental Systems Research Institute [ESRI], Redlands, California, USA), spatial analysis tools and images captured by the LandSat8 satellite operated by the United States Geological Survey (USGS; Reston, Virginia, USA). In each square, we chose a rice field at random, and in this field, we collected seeds using simple random sampling. Subsequently, seeds were stored at 4 °C for 3-mo.
Experimental design
The evaluation to the response of herbicides was carried out under controlled conditions, in a completely randomized design, with a complete factorial arrangement of three replicates. The first factor were the E. colona populations, along with the reference populations: susceptible (VM) and resistant (TA5). The second factor were the herbicides. Three herbicides were evaluated in three doses: bispyribac-sodium (sodium 2,6-bis((4,6-dimethoxypyrimidin-2-yl)oxy) benzoate; Army 400 SC, Adama Andina, Barranquilla, Colombia) at 0, 50 and 100 g ha-1; cyhalofop-butyl (butyl (2R)-2[4-(4-cyano-2-fluorophenoxy)phenoxy]propanoate; Clincher EC, Dow AgroSciences, Indianapolis, Indiana, USA) at 0, 270 and 540 g ha-1; and quinclorac (3,7-dichloroquinoline-8-carboxylic acid; Facet SC, BASF, Ludwigshafen, Germany) at 0, 375 and 750 g ha-1. The selected dose for each herbicide corresponded to 1x and 2x the commercial dosage. The experimental unit was a pot with 10 plants. To check out our results, the experiment was done twice.
The seeds were scarified with 10% HCl for 20 min, after that they were placed in Petri dishes with 0.6% agar (w/v), 0.2% KNO3 and 32 μg mL-1 ampicillin. When the seedlings emitted the first leaf, they were transplanted to soil in pots of dimensions 10×7×8.5 cm. Herbicides were applied when the plants displayed three fully extended leaves. Spraying was done using a CO2 at a constant pressure of 33 kPa with a Tee Jeet 8002 nozzle and an application volume of 200 L ha-1. Fourteen days after the treatment (dat) we evaluated fresh weight (g) and performed visual control (in contrast to the plants not treated with herbicide).
Statistical analysis
For the statistical analysis, variables were normalized using the following transformations: square root (bispyribac- sodium and cyhalofop-butyl) and logarithmic (quinclorac) for fresh weight, and arcsine for the control. Using the lme function of the nlme package with R Project (Pinheiro et al., 2017), we selected a linear mixed model using restricted maximum likelihood to analyze the three variables. We evaluated the inclusion of the populations as a random effect, as well as the interaction of the populations with the environment (experiments 1 and 2) and the replicates. The herbicide dose was included as a fixed effect in the model. In order to further evaluate the explanatory power of the model, we used a log-likelihood ratio test (Pinheiro et al., 2017), which is calculated as 2[log-likelihood of model B - log-likelihood of model A], where model A is more general than model B (model A has fewer random effects). The resulting distribution can be approximated as a chi-squared distribution with degrees of freedom KB-KA, where Ki is the number of parameters model i would predict (Pinheiro et al., 2017). We selected the simplest model, which resulted in significant results in the log-likelihood ratio test and the lowest values for the Akaike information criterion (AIC), and the Bayesian information criterion (BIC) (Aho et al., 2014). Subsequently, we compared the response of the variables (transformed) to the fixed effects, for different treatments, using the glht function of the multcomp package in R (Hothorn et al., 2008). Finally, for analyzing the results, average values obtained for each variable and each herbicide was using the original units.
Screening for resistance herbicides and spatial distribution
In order to categorize a population as herbicide resistant, we adapted the criterion of Panozzo et al. (2013; 2015b) in the following way:
where,
δrw
is relative fresh weight delta;
x¯wT
is average fresh weight of treated plants;
SDT
is standard deviation fresh weight of treated plants; is standard deviation fresh weight of untreated plants.
where
δrw
is visual control delta;
δC
is visual control delta; is average visual control of treated plants; SD is standard deviation visual control of treated plants.
S and SR categories were considered as susceptible and R and RR categories as resistant. When the categories were different within the populations according to the evaluated variable, the category obtained with the relative fresh weight was selected.
For each herbicide, the resistance categorization was represented spatially with the help of thematic cartography, using the ArcGIS 10.4.1 trial program, spatial analysis tools and images captured by the LandSat8 satellite operated by USGS. In order to establish if there was correlation between geographic location of the populations (field) and resistance, we computed spatial autocorrelation using a join-count analysis in the program R (R Core Team, 2014), assigning a value of 1 to resistant populations (R and RR) and 0 to susceptible populations (S and SR).
RESULTS AND DISCUSSION
The areas of the fields varied between 4 and 30 ha. In general, little information was obtained from herbicide applications. Of the information collected, 69% routinely used propanil, 46% bispyribac-sodium and 54% quinclorac. Specifically with the use of quinclorac and cyhalofop, it was established that their main use was in late post-emergence (beginning of flowering of E. colona) and in E. colona patches with low controls of previous applications. In addition, it was established that bispyribac-sodium is routinely used in each harvest cycle (2.3 cycles yr-1 in this area), with doses even 4 times higher than the commercial dose.
The mixed model that displayed the best fit, for each variable and each herbicide, is described in Table 1. It should be noted that general response to the herbicides (fixed effect) varied across populations (random effect) and was consistent among replicates (due to lack of interaction between populations and replicates). When the experiments were repeated after months, the response varied (random effect), but the variation was nonsignificant (due to the absence of interaction between fixed effects and environment; data not shown).
Bispyribac-sodium (Figure 1) did not diminish biomass accumulation (5%) and in fact, the response was similar even upon doubling the dose (p ≥ 0.05). The control increased by 12% upon doubling the dose (p ≤ 0.05). For cyhalofop-butyl and quinclorac the response varied with dosage (p ≤ 0.05). According to categorization of resistance (Table 2), 91% of populations were resistant (R and RR) to bispyribac-sodium, 48% to cyhalofop-butyl, and 43% to quinclorac (Figure 2). Additionally, we observed that 96% of populations became resistant to at least one herbicide (74% with 2x dosage), 65% to two (13% with 2x dosage) and 22% to three herbicides (0% with 2x dosage) (Table 2). The spatial autocorrelation analysis revealed that the resistance to all three herbicides had a random spatial distribution, and consequently there are not noticeable aggregation patterns on the maps (Figure 3).
Management problems with propanil-resistant E. colona in the 1980s led to the widespread use of bispyribac-sodium in Colombia (Valverde, 2007). In addition to the overuse of this herbicide encountered during sampling, in other research was established that bispyribac-sodium is the herbicide most used in the late post-emergency (Ramírez and Plaza, 2015). It is very likely that the strong selection pressure lead to the results we saw (91% of populations resistant). Additionally, and despite the fact that farmers in Colombia have a long history of difficulty managing E. colona with ALS-inhibitors, it was only in 2015 that Carranza and Plaza (2015) reported outbreaks of E. colona resistant to penoxsulam in rice fields of Tolima. These facts, along with the established history of using ALS-inhibitors in this area and our finding that the TA5 (resistant to penoxsulam) was also resistant to bispyribac-sodium suggest possible cases of crossed resistance, which still await verification.
Another alternative for the management of propanil-resistant E. colona were ACCase-inhibitors, mainly fenoxaprop, which with the subsequent selection of biotypes resistant to this herbicide, made way for an alternative herbicides such as cyhalofop-butyl (Valverde, 2007). In fact, bispyribac-sodium, propanil, and cyhalofop-butyl are herbicides frequently used in post-emergence in that area (Ramírez and Plaza, 2015). Cyhalofop-butyl was more effective at controlling weeds than bispyribac-sodium in the populations we studied (36% more for 1x, 49% more for 2x). Despite this, using the commercial dose, 48% of populations were resistant, doubling the dose only 9% of the populations remained resistant (RR). In other words, the problems in managing E. colona with cyhalofop-butyl, were overcome when the dosage was doubled. Both herbicides, ALS-inhibitor and ACCase-inhibitor, pose a high risk in the evolution of herbicide resistance (Valverde et al., 2000). In fact, applying a mixture of penoxsulam and cyhalofop-butyl, in E. crus-galli has been shown correlation with the selection of multiple resistance to these herbicides (Chen et al., 2016). The differences that we found when doubling the dose are due to historical selection pressure (changes in the use of cyhalofop-butyl when bispyribac- sodium entered the market; Valverde, 2007). Additionally, during the sampling process we found that unlike bispyribac- sodium, farmers do not use cyhalofop-butyl in all the crop cycles and when they do use it, they apply it in the late stages of post-emergence.
Table 1 Explanatory power and importance of the random effects in the mixed linear models selected for each variable with the bispyribac-sodium, cyhalofop-butyl and quinclorac herbicides.

x: Random effects adjusted in the selected model; P: populations; E: experiments; R: replicate; Loglik: logarithmic likelihood tests; AIC: Akaike information criterion; BIC: Bayesian information criterion; PRV: likelihood ratio test; df: degrees of freedom; p-value: significance between the models, based on Loglik.

A) Bispyribac-sodium at 0, 50 and 100 g ha-1, B) cyhalofop-butyl at 0, 270 and 540 g ha-1 and C) quinclorac at 0, 375 and 750 g ha-1. Variables: fresh weight and control. The average values are presented in the upper part of the bar, and in the lower part the significance of the comparison tests (different letters are significantly different, p ≤ 0.05).
Figure 1 Fresh weight and percentage of control of 23 populations of Echinochloa colona in response to bispyribac- sodium, cyhalofop-butyl and quinclorac.
The results with quinclorac were very similar to those of cyhalofop-butyl. This contrast with the difficulty that farmers report from the fields when they try to manage E. colona with quinclorac, whose application in the studied area is usually done in the last phase of post-emergence. Even though the product label does not specify the ideal moment to apply the herbicide during post-emergence, these results show that the use of the herbicide in the phase of three to four leaves could help diminish the difficulty of managing this weed species in the sampled fields.
Ten years after the introduction of quinclorac, O. Schmidt (2000 BASF, Germany), in personal communication with Valverde and Itoh (2001), reported quinclorac resistance in Colombia for the first time. The sole other study was that of Carranza and Plaza (2015), where 60% of the accessions were resistant. In our work, the percentage of quinclorac- resistant populations was 43%, but upon doubling the dose, only 9% were resistant (RR), meaning that similarly to cyhalofop-butyl, doubling the dose diminished the difficulty of managing E. colona. Even though these results suggest that increasing dose would be a good strategy to overcome the hardship in controlling E. colona, such an approach would also increase the selection pressure and would thus induce a more rapid evolution of the resistance to these herbicides. Cultural management to reduce the settlement, impact and fertility of these weeds in the crops should be the primary tool.
As far as the spatial distribution of the resistance is concerned, this was random for the three herbicides (there was no apparent pattern of aggregation). In other works on grass weeds resistant to ALS- and/or ACCase-inhibitor herbicides, which tried to establish some kind of correlation between geographic location and resistance, the results have been similar (Menchari et al., 2006; Schulz et al., 2014; Babineau et al., 2017). There have even been studies on the frequency of genes that confer resistance to ACCse in Alopecurus myosuroides, but they failed in detecting spatial structure in the distribution at all geographic scales (Menchari et al., 2006). Kuester et al. (2015) agree that resistance evolved independently in each field, meaning that control of the species is more successful if done at the level of each field.
Table 2 Categorization of the resistance of 23 populations from Echinochloa colona to bispyribac-sodium, cyhalofop-butyl and quinclorac.

R: > 1x (recommended dose) resistant; RR: > resistant to 1x and 2x (twice the recommended dose); S: > 95% susceptible with 1x; SR: > 80%-95% susceptible with 1x; light gray: > resistance to an herbicide; medium gray:
> resistance to two herbicides; dark gray: > resistance to three herbicides.
These results agree with those obtained through mixed model analysis, where introducing populations as random effects resulted in a response to herbicides that was differential across populations but homogeneous in each population. In other words, individuals in each field we sampled showed homogeneity in their response to treatments, but across fields, their responses were heterogeneous.
The selection pressure herbicides exercise locally, can explain the homogeneity of the populations at the level of each field. Due to, herbicide spraying programs are designed at the level of a field that create a mosaic of local selective pressures (Menchari et al., 2006; Délye et al., 2010). Additionally, when resistant genes are introduced into a field by pollen migration or resistant seeds from near-by fields (Busi et al., 2011), they most likely arrive repeatedly from different directions. Thus, spatial distribution of resistant plants in a field are expected to be patchy during a few generations from its introduction (Jasieniuk and Maxwell, 1994). Nybom (2004) found that for annual, autogamous plants that have gravity-dispersed seeds or tend to adhere to animals, a smaller genetic variance is expected intra-populationally, while inter-populationally this variance is likely higher. Echinochloa colona, being an annual and predominantly autogamous plant, is also awnless (Masood et al., 2016), which diminishes its mobility inside the field. In agreement with these aforementioned characteristics, the homogeneity among populations inside a field and heterogeneity across fields can be due to both selection pressure and to the biology of the species.

A) Bispyribac-sodium, B) cyhalofop-butyl and C) quinclorac, calculated from relative fresh weight and control. Categories: resistant (R and RR) and susceptible (S and SR).
Figure 2 Frequency of resistance of 23 populations of Echinochloa colona to, bispyribac-sodium, cyhalofop-butyl and quinclorac.
As far as the dispersion of the resistance is concerned, we must contemplate the ways in which resistant genes can propagate. In Costa Rica, Rojas and Agüero (1996) found that E. colona thrives in irrigation canals, dispersing seeds in the water and re-infecting commercial rice fields. The Saldaña and Purificación municipalities have been able to rely on irrigation water supplies from the Usosaldaña irrigation district during the last five decades, which could have favored the spreading of resistant seeds in the area. However, it is known that seeds can also be spread through the displacement of dirty agricultural machinery, the use of contaminated seeds, or by adhering to skin, feathers, animal paws and human clothing (Masood et al., 2016). During the sampling process, we observed that many agricultural activities were carried out by groups of persons and machinery operators that work in different fields, which would imply an anthropic dispersion of seeds in addition to that evidently achieved through irrigation.
On the other hand, we can establish that the timing of the experiments (3-mo gap between replicates) affected the response of the plants to the herbicides, as there was environmental variation between the two trials. The ability of groups of individuals present in the sampled populations to produce more than a single phenotype when exposed to environmental variation is consistent with the definition of phenotypic plasticity (Kelly et al., 2012). Nonetheless, the interaction between environment (experiments) and dose of herbicide (fixed effects) was nonsignificant; this means that the plastic response was not present in all the populations.
Starting from what we found with categorizing herbicide resistance and determining its spatial distribution, we argue that it is important to conduct any mitigation work at the level, of a single field or group of fields. Simple procedures such as cleaning machinery, clean crop seed, weed seed removal at harvest and seed destruction post-harvest, could help to decrease dissemination of seeds across fields (Busi et al., 2011), causing difficulty in the spreading of the resistance. Neve et al. (2009) signal the need to manage weeds using a broad perspective that takes into account the evolutionary biology of the plant species, achieving prevention or at least minimizing its adaptation and propagation. In this sense, these authors show the importance of a deep knowledge on the extent, structure and significance of the genetic variation inside and across weed populations. The use of mixed effects in this work allowed us to model variation among individuals, intra- and inter-populationally. Nonetheless, additional work is required in order to define this variation at a genetic level.

A) Bispyribac-sodium (recommended dose [1x] 50 g ha-1; twice recommended dose [2x] 100 g ha-1), B) cyhalofop-butyl (1x: 270 g ha-1; 2x: 540 g ha-1) and C) quinclorac (1x: 270 g ha-1; 2x: 540 g ha-1) (Scale 1:125 000).
Figure 3 Distribution map of the resistance of 23 populations of Echinochloa colona to bispyribac-sodium, cyhalofop- butyl and quinclorac.
CONCLUSIONS
Through this study, we documented herbicide resistance of Echinochloa colona in Colombian rice crops. We detected E. colona populations that were resistant to bispyribac-sodium (91%), cyhalofop-butyl (57%) and quinclorac (52%). In the case of cyhalofop-butyl and quinclorac, the percentage of resistant populations decreased upon doubling the dose of herbicide (13% and 17%, respectively). Ninety six percent of populations were resistant to at least one herbicide, and we detected multiple resistance to two (65%) and three (22%) herbicides. The spatial distribution of the resistance was random.
The inclusion of the mixed model in analysis enabled us to evaluate the response to the herbicides without allowing the natural variability of the populations to interfere with the results. This approach demonstrated that population responses were homogeneous inside the same field, heterogeneous across fields, and varied as a function of environment.
We conclude that the interaction between weed management and crop and biology of E. colona dictate the evolution of herbicide resistance at the level of a single field, and thus resistance management of this species of weed must be done at this scale too.