Immune and environment-driven gene expression during invasion: An eco-immunological application of RNA-Seq

Abstract

Host-pathogen associations change rapidly during a biological invasion and are predicted to impose strong selection on immune function. It has been proposed that the invader may experience an abrupt reduction in pathogen-mediated selection ("enemy release"), thereby favoring decreased investment into "costly" immune responses. Across plants and animals, there is mixed support for this prediction. Pathogens are not the only form of selection imposed on invaders; differences in abiotic environmental conditions between native and introduced ranges are also expected to drive rapid evolution. Here, we use RNA-Seq to assess the expression patterns of immune and environmentally associated genes in the cane toad (Rhinella marina) across its invasive Australian range. Transcripts encoding mediators of costly immune responses (inflammation, cytotoxicity) showed a curvilinear relationship with invasion history, with highest expression in toads from oldest and newest colonized areas. This pattern is surprising given theoretical expectations of density dynamics in invasive species and may be because density influences both intraspecific competition and parasite transmission, generating conflicting effects on the strength of immune responses. Alternatively, this expression pattern may be the result of other evolutionary forces, such as spatial sorting and genetic drift, working simultaneously with natural selection. Our findings do not support predictions about immune function based on the enemy release hypothesis and suggest instead that the effects of enemy release are difficult to isolate in wild populations, especially in the absence of information regarding parasite and pathogen infection. Additionally, expression patterns of genes underlying putatively environmentally associated traits are consistent with previous genetic studies, providing further support that Australian cane toads have adapted to novel abiotic challenges.

Keywords: Bufo marinus; cane toad; compositional data analysis; enemy release hypothesis; invasive species; spatial sorting.

Figure 1

Geographic distribution of the cane toad in Australia (dark gray region). Since arriving in Queensland in 1935, cane toads have further expanded their range through New South Wales, the Northern Territory, and into Western Australia. Black diamonds indicate our toad collection sites (from east to west): QLD (Gordonvale and Daintree, N = 5 each), NT (Cape Crawford and Timber Creek, N = 4 each), and WA (Caroline Pool and Durack River, N = 5 each). Map adapted from Tingley, et al. (2017)

Figure 2

Six unique patterns of gene expression in spleen tissue from invasive cane toads (Rhinella marina). We collected samples in populations from areas spanning the invaded range in Australia (QLD = Queensland, NT = Northern Territory, WA = Western Australia). Color indicates membership values of genes to clusters (purple = 0.7–0.8; pink = 0.8–0.9; red = 0.9–1). Tick marks on the x‐axis indicate sites across the toad's Australian range in which spleens were collected (Figure 1)

Figure 3

REVIGO plots displaying gene ontology (GO) terms (concepts/classes used to characterize gene function) depicting biological processes associated with transcripts following six major expression patterns in cane toads (Rhinella marina) across their Australian range (Figure 1). RNA‐Seq data from spleens were used to identify differentially expressed transcripts between invasion phases, then soft clustering was performed to visualize the expression patterns that these transcripts follow (Figure 2). Circles represent GO terms; those with the highest statistical significance are labeled. Circle size relates to breadth of GO terms. Colors show log₁₀ p‐values

Figure 4

Correlation between geographic distance and gene expression distance of invasive cane toad (Rhinella marina) populations across their Australian range (Figure 1). Euclidean distances in geographic location and gene expression between populations were calculated using the dist function in R. A mantel test (performed with the ade4 package) confirmed that these were significantly correlated (p = 0.005)