MATERIALS AND METHODS
2.1 Population sampling
We sampled populations of the mud-tidal snail Batillaria
attramentaria (G. B. Sowerby I, 1855) from two sites in Japan and one
site in the USA (details in Appendix 1: Table A1). We also include here
published data which we previously collected from a South Korean
population (Ho et al., 2019a). Population locations comprised Hacheon,
Cheollabuk-do, South Korea (on June 2016 at 35°32’N, 126°33’E, Ho et
al., 2019a); Nemuro city, Hokkaido Prefecture, Japan (June 2017 at
43°15’N, 145°28’E); Matsushima Bay, Miyagi Prefecture, Japan (May 2018
at 38°22’N, 14 1°4’E); and Monterey Bay, Elkhorn Slough, CA, USA
(February 2017 at 36°49’N, 121°45’W) (Figure 1). For each of the three
native populations (one Korean and two Japanese sites), 100 individuals
were collected. These native sites all had high surface salinities of
29-33 PSU at the time of collection. For the introduced (USA)
population, 50 individuals were collected. This site had a low surface
salinity of 4 PSU at the time of collection.
Due to river discharge, topography, and tides (e.g., Yoon & Woo, 2013),
different estuaries can experience wide fluctuation in daily salinity or
very little fluctuation. We characterized salinity fluctuation profiles
for each sampling site as low (Nemuro-Japan site; 27-34 PSU), moderate
(Hacheon-Korea site; 16-30 PSU), or high (Elkhorn Slough-USA site; 0-30
PSU) based on publicly available data on salinity fluctuation collected
over the past several years. These data were obtained from
http://www.nemuro.pref.hokkaido.lg.jp for Nemuro (Japan),
http://www.khoa.go.kr for Hacheon (Korea), and http://www.mbari.org for
the Elkhorn Slough (USA). We could not find recent records of salinity
fluctuation at Matsushima Bay (Japan); however, past data on average
monthly salinity levels inside Matsushima Bay indicates that salinity at
this site fluctuates from 27 to 34 PSU (Ventilla, 1984) and is similar
to the Nemuro site.
In animal locomotion experiments, it can sometimes be advantageous to
use individuals of similar size. However, we noticed in the field that
the typical body size of the introduced population was larger than the
native populations. Since our primary question centers on how different
populations respond to salinity stress, and since body size might play a
role in both locomotion and salt tolerance, we allowed body size to
differ between the native and introduced populations. This decision
somewhat constrains our ability to conclude whether body size is a
driving factor behind salinity responses; however, it allows us to use a
representative sample of each population rather than using specimens
whose size is not typical of their population, and thus makes our study
more ecologically relevant. Ultimately, the native and introduced
samples we collected differed by about 1 cm (average native shell length
= 2.1 cm; introduced shell length = 3.1 cm). All collected specimens
were maintained in a plastic aquarium with a constant air temperature of
25°C, a water salinity that was the same as their collection site, and a
12h L:12h D photocycle for two days prior to the experiments, in order
to reduce the effects of transportation stress and to allow for
acclimation to the lab.
2.2 Salinity stress
experiment
For each population consecutively, we conducted a 30-day experiment
examining locomotor behavior under different salinity treatments. We
randomly divided each population into five treatment groups, with 20
individuals per group for the native populations and 10 individuals per
group for the USA population. These groups were maintained in separate
plastic aquaria (40 × 23 × 21 cm3, with an inclined layer of sea sand
set up on the bottom and one liter of aerated artificial saline water,
Supplemental Figure S1, Ho et al., 2019b) at salinities of 43, 33, 23,
13, and 3 PSU for 30 days. Saline water was freshly prepared every two
days from overnight-aerated distilled water and Instant Ocean Sea Salt
(United Pet Group Inc., Cincinnati, OH, USA). Snails in each group were
marked with nail polish (Eco Nail color, Innisfree, South Korea) to keep
track of their identity. All animals were fed to satiation every two
days with a commercial brand of fresh, chopped seaweed (Ottogi, South
Korea) throughout the 30-day experiment.
2.3 Locomotor behavior tracking
We recorded the movement of each snail for one hour every two days
throughout the 30-day experiment following the protocol in (Ho et al.,
2019a). Briefly, we used a Sony NXCAM camera (AVCHD Progressive MPEG2
SD) to film individual snails in the center of a disposable Petri dish
(diameter: 9 cm) filled with artificial seawater which had the same
salinity as the snail’s assigned treatment group. We increased the video
playback rate using AVS Video Editor v.7.1.2.262 and cropped
videos using Avidemux v.2.6.12. We used the Spectral
Time-Lapse (STL) toolbox (Madan & Spetch, 2014) implemented in
Matlab release R2014a (MathWorks Inc.,
Natick, MA, USA) to quantify movement distance of the snails.
2.4 Shell measurements
We measured shell length of all individuals using images extracted from
the recorded videos. We used K-Multimedia Player software (KMP Player,
PandoraTV, Pankyo, Korea) to extract one video frame (resolution about
1200 × 1200 pixels) for each snail, and then used tpsDig2 to digitize
the most anterior and posterior points of the shell (Rohlf, 2008;
Appendix 1: Figure A1). The diameter of the petri dish (9 cm) was used
as the conversion scale to calculate snail length.
2.5 Sequencing and phylogenetic
analysis
After the salinity experiments, we extracted genomic DNA from fresh foot
tissue of all snails using a Dneasy Blood & Tissue kit (Qiagen, Hilden,
Germany). We PCR-amplified the mitochondrial CO1 gene using
published CO1 primers (Ho et al., 2015) and a Fastmix/Frenchetm
PCR kit (IntronBio, South Korea). PCR products were purified using a Dr.
Prep kit (MGmed, Seoul, South Korea). Sequencing reactions were
performed using a Bigdye Terminator V3.1 Cycle Sequencing kit (Bionics;
Seoul, South Korea).
We performed a Bayesian phylogenetic analysis based on the CO1sequences of five representative specimens which were identified as
Kuroshio and Tsushima haplotypes and from Korea, Japan, and the USA
(GenBank accession no.: MG241503-06 and MT800763), and 53 previously
sequenced specimens from the shorelines of Korea (GenBank Accession: No.
HQ709362– 81, Ho et al., 2015) and Japan (GenBank Accession: No.
AB164326-58, Kojima et al., 2004). We used a closely related species,Batillaria multiformis, as an out-group (GenBank Accession: No.
AB054364, Kojima, Ota, Mori, Kurozumi, & Furota, 2001). We employed
MegaX (Kumar, Stecher, Li, Knyaz, & Tamura, 2018) to predict
the best substitution model for the CO1 data, and found that the
Hasegawa-Kishino-Yano + Gamma + Invariable (HKY85+G+I) substitution
model was the best fit based on its corrected Akaike information
criterion (AICc) value of 1500.65. We then ran the analyses applying
Maximum Likelihood statistical method (Nei & Kumar, 2000) using
MegaX with 1000 replications of bootstrap and Neighbor-Joining
statistical method (Saitou & Nei, 1987; Studier & Keppler, 1988) using
Geneious tree builder incorporated in Geneious
v.6.1.8 (Kearse et al., 2012) with 1,000,000 of random seed to
construct a phylogenetic tree.
2.6 Model terms
Our main purpose was to assess the effect of differential salinity
exposure (es) on the locomotion of B. attramentaria (see
the description of the 30-day experiment above). In addition to
es, we also included in our analyses four other model terms
that might influence snail response to salinity: origin (o),
location (lo), population (p), and CO1 lineage
(li). Origin was defined as either native (pooled three
populations from Korea and Japan) or introduced (one population from the
USA). The location term indicated the countries where the snails were
collected (Korea, Japan, or the USA). The population term described the
sampling site (Hacheon in Korea, Nemuro city and Matsushima bay in
Japan, and Elkhorn Slough in the USA). Lineage was defined as either
Tsushima (comprising the Hacheon and Nemuro populations) or Kuroshio
(comprising the Matsushima and Elkhorn Slough populations) based on the
individual’s position in the CO1 phylogenetic tree.
2.7 Statistical analyses
For all four populations, snails exposed to 3 PSU did not move at all
and died within the first 16 days of the experiment, apparently due to
the extreme osmotic stress represented by such low salinity. Locomotor
data for all 3-PSU groups were therefore excluded from all analyses. All
other individuals survived the entire duration of the experiment and
were included in the analyses.
We applied a linear mixed-effect model (LMM) to assess the impacts of
multiple predictor variables (see “Model terms” above) on the
locomotor response of B. attramentaria to different levels of
salinity stress using the package nlme version 3.1-140 (Pinheiro, Bates,
DebRoy, Sarkar, & R, 2018) implemented in R version 3.0.2 (R
Development Core Team, 2011). Since we wished to investigate
specifically the impacts of salinity stress, geographic distribution,
and genetic composition, we fit a set of competing models with separate
single predictor variables (es, o, t, lo, p, and li)
and their interactions (es × o, es × o + li, es × lo, es × lo +
li, es × p, es × p + li, es × li, es × li + o, es × li + lo, and
es × li + p). We treated snail identity as a random variable.
For the sake of simplicity, we limited the maximum number of predictor
variables to three per model. Since body size closely corresponded to
the origin factor, which was already in our LMM analysis (i.e., we knew
native snails were smaller than introduced snails), we did not include
body size as a factor in the LMM analysis. We centered all the predictor
variables to mean = 0 and standardized to s.d. = 0.5 to remove
multicollinearity and to directly interpret the results in terms of
effect size (allowing us to compare predictors). The response variable
was the movement distance of snails, which was measured every two days
throughout the 30-day acclimation experiments. Prior to the analysis, we
square-root transformed the movement data using the package rcompanion
2.2.1 (Mangiafico, 2017) implemented in R 3.0.2 (R Development Core
Team, 2011) to meet the assumption of normal distribution. To determine
the best covariance structure for the LMM tests, we tested our response
variable against several covariance structures: First Order
Autoregressive (AR1), Compound Symmetry (CS), and Unstructured (UN). We
compared their corrected AICc values using the package MuMIn 1.9.13
(Barton, 2009) for R 3.0.2. The AR1 model was the best-supported
covariance structure based on the AICc value (3407.24 versus 3710.22
(CS) and 3708.22 (UN)) and was therefore chosen for the LMM tests as the
best available compromise between bias and lack of precision.
We next applied a multimodel inference procedure (Burnham & Anderson,
2003) to the set of competing linear mixed-effect models (LMMs) to
select the most parsimonious model that best described snail locomotor
response. The models were compared based on AICc values using the
aforementioned MuMIn package. The model with the lowest AICc value and
those satisfying a ΔAICc ≤ 6 cut-off rule (Richards, 2005, 2008) were
considered the best-fit or most parsimonious models. We then conducted
post-hoc multiple comparison tests of the best-fit model to examine the
effects of each explanatory factor. Additionally, we used MuMIn to
perform model averaging and estimate the importance of predictor
variables by summing the weights of models where the variables appeared.
In parallel, we examined differences in shell length between native and
introduced groups and among native populations using a two-way ANOVA.
The first level involved comparisons of snails from four populations,
three locations (i.e., countries), and two origins (introduced vs.
native). The second level involved comparing individuals belonging to
the two mitochondrial CO1 lineages of Kuroshio and Tsushima.
2.8 Data deposition
Data available from the Dryad Digital Repository:https://datadryad.org/stash/dataset/doi:10.5061/dryad.455mv2mand the Mendeley Data repository:http://dx.doi.org/10.17632/jjjmh26c2g.1.