Seasonal plasticity in sympatric <em>Bicyclus</em> butterflies in a tropical forest where temperature does not predict rainfall

While temperature variation appears to be the main environmental cue for plasticity in adult traits in many species of Mycalesina, relying on temperature would result in a mismatch between adult phenotype and environment in some regions. We measured phenotypes of six species of Bicyclus butterflies (Nymphalidae: Satyrinae: Mycalesina) in a humid tropical forest with two rainy seasons per year and modest unimodal seasonal temperature variation, such that temperature does not predict rainfall and butterflies can reproduce year-round. The butterflies showed subtle temporal variation in body size and relative eyespot size, while relative androconia length was robust to temporal environmental variation. After higher temperatures, body size tended to be smaller, and relative eyespot size was larger for some species-eyespot combinations. This indicates that these butterflies follow the “hotter is smaller” rule and show developmental plasticity in eyespot size that is typical in this clade. Eyespot sizes tended to be correlated with each other, except Cu1 in B. auricruda and some eyespots that always remained very small. Androconia length was not related to eyespot size. This pattern of correlations suggests conserved cue-use and shared mechanisms for eyespot size using both temperature and rainfall-related cues, with some exceptions. Methods Data collection Our study site was a sub-montane tropical forest near the Makerere University Biological Field Station (0°13’ - 0°41’N and 30°19’ - 30°32’E) in Kibale National Park, Western Uganda. In this region, there are two rainy seasons per year, while temperature has a unimodal distribution so that there is a warm and a cool dry season (Valtonen et al. 2013). Thus, unlike in study sites of previous studies (Oostra et al. 2014a, van Bergen et al. 2017), there is probably selection against developing a dry-season phenotype when temperatures are lower (and vice versa), at least during part of the year. During the study period, a data logger (Lascar EL-USB-2-LCD) was placed inside the forest, suspended 1 meter above the ground, and protected from direct rain and sunlight by a plate of roofing zinc. The logger recorded temperature and relative humidity in half-hour intervals. In addition, daily rainfall and temperature data were obtained from a weather station at the field station (within 1km of all butterfly sampling locations; Chapman et al. 2018). To minimize the impact of our study on the local butterfly populations, we focused our sampling efforts on male specimens. Up to ten individuals of the six most abundant Bicyclus species (B. collinsi (Hewitson, 1873); B. mollitia (Karsch, 1895); B. smithi (Aurivillius, 1899); B. auricruda (Butler, 1868); B. golo (Aurivillius, 1893), and B. graueri (Rebel, 1914)) were collected weekly from baited traps for 14 consecutive months (23rd July 2013 to the 26th of September 2014). All six species have eyespots in the distal region of the ventral wings (Figure 1). Whilst some species have multiple androconial brushes, all the investigated species have a prominent brush with its base located in the dorsal wing cell of the hindwing. This shared brush was selected as the androconial trait that could be compared across species. All four wings of collected butterflies were placed on a Nikon grey card, which was placed on graph paper (Figure 1) and photographed using a Nikon D7000 camera in a custom-made studio with constant light conditions (luminance and intensity) and the same manual settings of 1/125 shutter speed and F14 aperture for all photographs. We then used a macro in ImageJ to measure a proxy of the wing area of each wing, the area covered by four ventral eyespots, and the length of the basal hindwing androconial brush (Figure 2). Data analysis To obtain the average temperature for each day, we averaged the maximum and minimum temperatures provided by Chapman et al. (2018) rather than our own measurements, as these also cover the months before butterflies were collected (necessary for cross-correlation analyses). We averaged rainfall and humidity by two-week period and month. As a proxy of body size, we averaged the area of triangles measured from forewings and hindwings (Figure 2). In general, wing area is a well-established proxy for body size in these butterflies (e.g. Bergen et al. 2024). We calculated relative eyespot size as eyespot area divided by the body size proxy. Relative androconia length was calculated as the length of androconia divided by wing length. To avoid periods with missing data due to low abundance of butterflies, species’ averages of traits were calculated per two weeks or four weeks, depending on species abundance. The few remaining missing data points (5 out of 135) were replaced by the average of data points from two weeks before and two weeks after (for biweekly data) or four weeks before and four weeks after (for four-weekly data). To assess whether butterflies showed seasonal dimorphism (distinct wet season and dry season morphs), we generated density plots of body size, eyespot size, and androconia length for each species using the R function geom_density from the package ggplot2 (Wickham 2016, R_Core_Team 2024). We visualized temporal trends using the loess function in R (R_Core_Team 2024). To estimate whether traits were linked to habitat seasonality (seasonal changes in average phenotype), we performed autocorrelation analyses using the function ggAcf in the R package ggplot2 (Hyndman & Khandakar 2008). Autocorrelation analysis tests whether, within a single time series, there are correlations between data points that are a particular time lag apart, as would be the case with seasonal patterns (e.g, correlation between data points that are six months apart). To focus on seasonality, for each autocorrelation, we visually identified the lag with the most negative correlation coefficient and determined whether it was statistically significant, thus ignoring lags next to zero that typically show positive correlation coefficients. To test if temperature and rainfall drive temporal variation in butterfly traits, we performed cross-correlation analyses using the function ggCcf in the R package ggplot2 (Wickham, 2016). In addition, we performed such analyses for relative humidity (related to both rainfall and temperature). Cross-correlation analysis tests whether two time series are correlated with each other with a certain lag. We expected the lag to range within the length of development time, i.e., time between egg hatching and adult eclosion (about 6 weeks; Molleman et al*.* 2016, van Bergen et al. 2017) and development time with an added month to account for adult life span (Molleman et al. 2007) and an extra month if the effect is mediated by host-plant growth (Valtonen et al., 2013; total 4 months lag). For each cross-correlation, we visually identified the lag with the highest correlation coefficient and noted the sign and whether it was statistically significant, focusing on lags of less than six months. Among butterfly species, we tested for correlation within traits (no lags), which may indicate an evolutionarily conserved mechanism across species. Within species, we tested for the degree of correlation among traits, which would indicate a linked developmental mechanism. Given the low number of species included in our sampling, we did not account for phylogenetic non-independence among species in our cross-correlation analyses.