Ben Sulser 0000-0002-8750-0942
· sulserrb
Evolutionary Ecology Group, Institute of Ecology and Evolution, University of Bern, Baltzerstrasse 6, 3012 Bern, Switzerland
· Funded by Bern Burgergemeinde
Sheila Christen
Catherine L. Peichel 0000-0002-7731-8944
· cpeichel
Division of Evolutionary Ecology, Institute of Ecology and Evolution, University of Bern, Baltzerstrasse 6, 3012 Bern, Switzerland
· Funded by Swiss National Science Foundation (TMAG-3_209309/1)
Ruslan Hlushchuk✉ 0000-0001-2345-6789
· RuslanHlushchuk
microCT research group, Institute of Anatomy, University of Bern, Baltzerstrasse 2, 3012 Bern, Switzerland
✉ — Correspondence possible via GitHub Issues
or email to
Ruslan Hlushchuk <ruslan.hlushchuk@unibe.ch>.
Abstract
Can we predict evolution?
The three-spined stickleback (Gasterosteus aculeatus) is a well-recognized system for understanding adaptation to divergent habitats.
Populations of benthic (benthos-feeding) and limnetic (water-column feeding) stickleback differ in a number of phenotypic traits that are associated with shifts in dietary specialization.
Modern investigation of the evolutionary changes in this organism often requires the analysis of hundreds, if not thousands, of individuals.
However, analyses of the structures required for feeding – especially the jaws and complex internal branchial anatomy — require considerable time and expertise, with destructive sampling and fine dissection skills needed for quantitative analysis.
The advent of X-ray microtomography and 3D-scanning technology affords non-destructive sampling and an increases the amount of high-resolution data available for study, but at the substantial cost of increasing complexity and processing time for each specimen such that these techniques are often unfeasible for studies of eco-evolutionary scale.
To address these concerns, we developed a rapid and semi-automated segmentation and analysis pipeline based on both the Jupyter interactive development environment and the Biomedisa image segmentation platform for investigating three-dimensional morphological adaptation within the three-spined stickleback.
The pipeline includes splitting multi-specimen scans into regions of interest for each specimen, reconstruction of targeted anatomy, and morphometric analyses.
We then applied this pipeline to a sampling effort comprising 44 multi-specimen scans representing 216 specimens of divergent benthic and limnetic stickleback populations, showcasing the possibility of using high-throughput scanning data to provide tests of ecological and evolutionary hypotheses.
Introduction
The three-spined stickleback (Gasterosteus aculeatus) is an oft-studied organism for understanding the independent evolution of similar traits in similar environments (Bell & Foster, 1994; Reid et al., 2021).
This species exhibits marked differences in marine-freshwater, lake-stream, and benthic-limnetic ecotypes (WILLACKER et al., 2010).
This study will focus on the benthic-limnetic axis, using samples from a long-term evolutionary experiment currently investigating divergent populations of limnetic and benthic stickleback within the Kenai Peninsula of Alaska (USA) (Hendry et al., 2024).
This project, the Forward In Time Natural Experimental Study of Selection (FITNESS), aims to study the predictability and repeatability of evolution.
Two pools of sticklebacks — one made from four source populations of limnetic and four source populations of benthic sticklebacks — have been placed into eight destination lakes, four of which are small and benthic and four of which are large and limnetic.
These new populations have been sampled every year in order to track the genotypic and phenotypic trajectories of these introduced populations.
Understanding the initial variation in the source populations is essential to this project, as this initial variation would be expected to reflect which phenotypes are associated with each ecotype under study.
Among these and other bony fish, differences in jaw structures are directly related to functional and kinematic differences between different ecotypes (Haines et al., 2020).
Benthic stickleback have modified jaws for enhanced suction force and hypertrophied epaxial muscles to aid in foraging on benthic invertebrates and by contrast, limnetic stickleback have modifications for larger jaw protrusions and quick strikes during ram feeding (McGee et al., 2013).
The internal hyoid arch-branchial arch complex is an important structure implicated in diet and feeding ecology (BERNER et al., 2008; Schluter & McPhail, 1992).
While the shape and arrangement of paired ceratobranchial and pharyngobranchial bones within this complex aid food processing and water vortex generation during feeding (Brooks et al., 2018), the shape of these bones have received comparatively little attention relative to other aspects of dietary anatomy.
This is likely due to the flattening and destructive sampling used in traditional raker counting methods, which dissect and deform these structures to render them visible for manual measurement (Ellis & Miller, 2016).
These structures are, however, difficult to study without full cranial dissection and corresponding distortion of the branchial anatomy.
3D analyses preserve these features at a high resolution.
This work is embedded within the Genomics axis of the Alaska Stickleback Restoration Project, with which Katie Peichel, Ben Sulser and Sheila Christen are affiliated.
Micro-computed tomography
X-ray microtomography (μCT) is an indispensable tool to gain non-destructive insights into the inner structure of highly diverse samples, namely for specimens studied in the biomedical sciences (Rawson et al., 2020).
Microtomographic imaging is ideally suited to non-destructively assess the morphology of different fish species (Ford et al., 2023), including the internal anatomy and small structures difficult to quantify without additional preparation.
While these structures can be rendered by hand by a skilled investigator, the time and cost required per-specimen is inefficient for the scale required via eco-evolutionary study and this requires destruction of the mandibular and cranial anatomy of the specimen.
This project aims to address these gaps, demonstrating a novel pipeline for rendering and auto-splitting of a multi-specimen scan for mass sampling, creating a dataset with consistent parameters that can be fed to downstream machine learning techniques (Lösel et al., 2020) to aid in the segmentation of individual bony structures in each scan.
Once a Biomedisa model is trained, the entire pipeline runs from multi-specimen input to rendered structures for each specimen in a fraction of the time and resources used in traditional analysis.
Materials & Methods
Figure 1: Workflow overview
Sample procurement and preparation
The specimens used for this study were collected from source lakes as a part of the FITNESS project in the region of Cook Inlet, Alaska.
Fish were collected using unbaited minnow traps in two separate field seasons, the first taking place from May 26–June 10 2023 and the second taking place from May 25–June 11 2024.
Specimen collections were taken from a random sample of 30 fish from each lake, under Alaska Department of Fish and Game (ADFG) permits SF2023-030 and P-24-015 for 2023 and 2024, respectively.
Fish were euthanized with MS-222, photographed, and preserved in 10% formalin in a bag with a specific label, under Animal Use Protocol (AUP) MCGL-8265.
At the end of each field season, samples were shipped from Anchorage (AK, USA) to Bern (BE, CH) where they were stored until scanning time.
Due to their inherent contrast difference to the surrounding tissue, the structures of interest in this study (teeth and bones, i.e., jaws and skull) are well visualized in unstained samples, hence no further preparation of the fish was necessary.
μCT imaging
In a small pilot study, we determined the optimal scanning parameters to meet the constraints on total scanning time, resolution and sample handling.
To optimize for these constraints, we scanned all the sticklebacks in batches of 6 fish in a custom-made 3D printed sample holder in a single scan.
This holder was generated with OpenSCAD and is available online, either directly as STL file for printing or as (parameterized) OpenSCAD file for adaptation to other classes of samples.
Both files are part of a library of 3D-printable sample holders for tomographic imaging (Haberthür, 2019).
Tomographic imaging was performed on a Bruker SkyScan 2214 at the Institute of Anatomy, University of Bern, Switzerland.
In total we performed 44 scans, each of the scan usually containing 6 fish in the sample holder.
The relevant details of each scan are collated in a table in the Supplementary Materials; a short overview of the scanning parameters is given below.
The X-ray source was set to a voltage of 60 kV and a current of around 110 µA for all but one scan where we used a source voltage of 49 kV and 159 µA due to operator error.
For each sample, we recorded a set of 3601 projections of approximately 3000 x 2000 pixels at every 0.1° over a 360° sample rotation.
Every single projection was exposed for about a second (depending on the sample).
Due to the length of the fish, we had to acquire so-called stacked scans, on average we scanned 3 fields of view along the rotation axis of the sample holder.
This resulted in scan times between 3 to 5 hours.
The projection images were then subsequently reconstructed into stacks of 8bit PNG images with NRecon (Bruker microCT, Kontich Belgium, Version: 2.1.0.1 or 2.2.0.6), without applying any ring artefact or beam hardening correction.
The isometric voxel sizes in the resulting datasets vary from 15 to 19 µm.
Data analysis
Preparation and handling of tomographic datasets
After acquisition, a simple script was used to copy the relevant data to both archival storage and storage accessible by all co-authors at the same time.
The preview notebook is used for surfacing issues with the scanning.
For this, we read all relevant scanning and reconstruction parameters from the log files of each scan.
Afterwards, we efficiently loaded the reconstruction PNG images from disk with the dask_image.imread.imread function (Dask Development Team, 2016).
Like so, we can map all the generated reconstructions to memory and quickly generate maximum intensity projections (MIP) of each scan (see Figure 2 for an example) for both quality control and further processing.
Figure 2: Maximum intensity projections of one acquired dataset along all three cardinal axes.
Separation notebook
The separation notebook processes all the performed scans to extract each individual fish from each scan encompassing 6 fish in total.
As in the preview notebook, we efficiently load all the PNGs from disk with dask(Dask Development Team, 2016).
Based on the previously extracted MIP images and a simple labeling of these images (skimage.measure.label), we extract both the labels in the custom-made sample holder and the positions of single fish in the scan (skimage.measure.regionprops) (see Figure 3).
This extraction is completely reproducible and well-adapted to the custom-made sample holder.
Figure 3: Automatically detected regions based on maximum intensity projection along the rotation axis of the tomographic scan.
The regions are numbered consecutively from the top left to the bottom right.
These numbers are mapped to the correct fish ID in the next step.
Based on a simple mapping of the detected region to the ID numbers of the scanned fish, we labeled the resulting images and presented these images together with photos of the lab book and sample tubes for verification (see Figure 4).
Figure 4: Mapping lab book notes, photos and detected regions to fish ID.
The skimage.measure.regionprops function we used for labeling not only returns the positions of each detected fish, but also the extent of the bounding box of each detected region.
We extracted each region of each fish separately out of the large reconstructions (with a configurable border buffer, see Figure 5) and wrote these extracted regions to disk in discrete folders for efficient further analysis.
In a first step, we wrote the regions of the single fish to disk in zarr(Miles et al., 2020) format, which is a preferred format to store n-dimensional arrays on disk.
In addition to this, we also wrote a log file for each extracted region, containing all relevant information to redo the cropping step completely manually (an example of such a log file is shown as part of the processing repository).
Figure 5: Double-checking crop extent and fish ID.
The top row shows the extracted regions from the previously calculated MIP of the full scan, the bottom row shows the MIP images of the extracted regions.
Both rows must show exactly the same region.
Writing the regions as zarr files made it possible to efficiently work with the image data of each extracted fish and to convert that data to any desired format for further analysis.
For this further analysis, we wrote stacks of PNG images and additionally, as nrrd files for each fish region as cropped as well as cropped and binarized regions out of the original dataset.
These binarized regions were segmented into bone and background based on a simple multi-level Otsu thresholding method (廖炳松 et al., 2001).
Figure 6: Three-dimensional preview of extracted region, automatically thresholded.
Extraction of features of interest
After separation, the cropped image files were checked and rendered via the use of 3D Slicer (Kikinis et al., 2013) and the SlicerMorph extension (Rolfe et al., 2021).
The individual elements of the branchial apparatus were rendered using a combination of thresholding and split islands tools to separate the pharyngobranchials, epibranchials, basibranchials, hypobranchials and ceratobranchials (see Figure 7).
Figure 7: Example of Branchial Anatomy with CB1 and CB2 highlighted. Abbreviations: PB = pharyngobranchials, EB = epibranchials, BB = basibranchials, HB = hypobranchials, CB = ceratobranchials.
Once rendered, these bones were exported as a colored labelmap alongside the .nrrd file from which they were segmented to pass to the Biomedisa program.
Machine learning and model training
As a group, a dataset of 51 specimens (including .nrrd and .label files) were passed to Biomedisa (Lösel et al., 2020) to train a segmentation model.
We allowed the for rotation of 180° to account for possible specimen variability, and an 80/20 split between training and validation data.
The model was trained with a batch size of 24 and 50 epochs, under Network architecture 32-64-128-256-512.
Landmarking of models
To demonstrate the effectiveness of this tool and the importance for 3D morphometrics for answering eco-evolutionary questions, we have run a demonstration quantifying the shape differences of the ceratobranchial bones.
Once trained, we applied the Biomedisa segmentation model to the remaining 160 specimen volumes and landmarked the final results using Stratovan Checkpoint (Stratovan Corporation, n.d.).
As a test and for subsequent analysis, the first two ceratobranchials on the right side were chosen for comparison across all specimens.
Type II landmarks were set on the ends of each bone, with semilandmarks in-between each to cover axes of curvature along the bone (see Figure 8).
In total 7 landmarks and 4 semilandmark curves (two containing 20 semilandmarks, two containing 15) on the first ceratobranchial (CB1), and 5 landmarks and 3 semilandmark curves (one containing 20 semilandmarks, two containing 15) on the second ceratobranchial (CB2).
Equal distances were ensured using the resample_curves function in 3D Slicer.
Figure 8: Landmarking shown on dorsal (a) and ventral surfaces (b).
Analysis of shape
All subsequent analyses were run using R (version 4.4.1, (R Core Team, 2021) and the geomorph package (Adams & Otárola‐Castillo, 2013).
Both bones were split and analyzed separately after generalized Procrustes analysis (GPA) using the gpagen() function, with Principal Component Analysis (PCA) and linear models run with gm.prcomp() and procD.lm(), respectively.
Linear fits were further investigated via the pairwise() function to analyze differences in pairwise statistics.
Results
μCT data
Acquisition and reconstruction of fish proved successful and efficient.
216 unique specimens were scanned in a total scanning duration of 18 days, 12 hours, and 6 minutes.
We acquired 158444 projections, reconstructed into a total of 177749 reconstructions, to about 4040 files per scan (N=44).
The total size of this sampling effort is ~44 GB of .zarr files, ~64 GB of .nrrd files.
Fish separation
Our method reproducibly extracts each of the 6 fish scanned simultaneously in one scan.
The custom-made sample holder aligns the single fish in the vertical axis around the rotation axis of the tomographic scan.
The extraction based on the MIP image along the rotation axis is completely automated and very robust, since the detected fish ‘regions’ do not overlap in the resulting image.
Depending on the available hardware, it may not even be possible to load the full stack of each scan into a software to manually perform the cropping, such as Fiji (Schindelin et al., 2012).
Large stacks of images (in other words larger than the RAM of the available machine) can be loaded as virtual stacks, but to manually crop the region of each fish from the large scan with the Crop (3D) function, one needs to load the full dataset.
Since one (exemplary) dataset (Sticklebucket_10) is 7 GB on disk and reported as being 35.4 GB when loaded in Fiji, using the 3D cropping function on an uncropped single dataset is not possible without using a powerful workstation.
Extracting the single fish from the encompassing dataset would thus be a two-step manual process, e.g. cropping the full dataset loaded as virtual stack and then cropping it down further before writing out the cropped stack.
For each encompassing scan this would need to be repeated 6 times (for each of the 6 fish in each of the encompassing scans).
In addition, such a manual process is not reproducible in the sense that it cannot be consistently replicated by others using the same data since the manual cropping is operator-dependent.
Algorithmically/automatically cropping the large datasets based on the axial MIP image leads to both reproducible cropped regions and efficiently uses the operator time (namely no operator time) (see Table 1).
Table 1: Estimates of time comparisons between manual and pipeline runs.
Task
Est. Manual Time [min]
Pipeline Time [min]
Scanning single scan
#
##
Splitting and rendering volumes
#
##
Segmentation
10-15
2-3
Our automated extraction process also writes human-readable log files documenting the cropping position in the encompassing dataset and the crop extent.
This enables reproducible double-checking and confirmation of the process after the fact (see this direct link for one such log file).
Thresholding
The separated fish were segmented based on a simple multi-level Otsu thresholding method.
This relatively simple segmentation was sufficient to extract all the features we analyzed further, and we did not have to employ more advanced thresholding methods in our separation pipeline.
Selection and individual rendering of the branchial structures takes between 10-15 minutes; the average Biomedisa render takes 2.5 minutes once trained (See Table 1).
Analysis
The speed and quality of these data allow us to study the internal branchial anatomy at scale and in situ, without the need for fine dissection.
While the shape and arrangement of the ceratobranchials and the corresponding bony gill rakers are hypothesized to work in tandem for food processing and water vortex generation during suspension feeding (Brooks et al., 2018), the shape of these bones has received comparatively little attention.
This is likely due to the flattening and destructive sampling used in traditional raker counting methods, which dissect and deform these structures to render them visible for manual measurement.
3D analyses preserve these features at a high resolution.
After GPA alignment, we are able to quantify the shape differences among all fish scanned for this project.
Changes due to allometry (using the metric of centroid size or standard length of the fish) were significant, but slight: explaining only a small fraction of shape variation in both bones.
Both linear models and PCA results suggest that the lakes themselves - and not overarching categories of ecotype or sex - drive most of the shape variation in these bones (CB1: p = .001, Rsq = 0.03246, CB2: p =.001, Rsq = .06220).
The ecological variation present across the first ceratobranchial shows a significant but quite small effect with lake origin (p = .009, Rsq = 0.02056 ), and with a large amount of overlap in the resulting shape space (see Figure 9).
Figure 9: PCA of CB1 colored by lake ecotype.
Warps are indicated at the extremes of each axis by vectors drawn from a mean shape (red).
The second ceratobranchial bone, on the other hand, shows equally small yet significant shifts associated with the ecotype (p = .001, Rsq = 0.0377).
The pattern of difference between benthic and limnetic gill rakers are, for this bone, clearly divergent in shape space (see Figure 10).
Figure 10: PCA of CB2 colored by lake ecotype.
Warps are indicated at the extremes of each axis by vectors drawn from a mean shape (red).
These differences in ecological patterning were also broken down by lake (See Figures 11 and 12).
Figure 11: PCA of CB1 colored by individual lake.
Warps are indicated at the extremes of each axis by vectors drawn from a mean shape (red).
Figure 12: PCA of CB3 colored by individual lake.
Warps are indicated at the extremes of each axis by vectors drawn from a mean shape (red).
The differences in the 2nd ceratobranchial appear to be driven by divergence in the South Rolly population, supported by significant pairwise differences observed between this lake and all other lakes observed in CB2 and not in CB1 (see supplementary information).
After GPA alignment, we are able to quantify the shape differences among all fish scanned for this project.
Changes due to allometry (using the metric of centroid size or standard length of the fish) were significant, but slight: explaining only a small fraction of shape variation in both bones.
Discussion
Pipeline and efficiency
Once all elements of the pipeline are together, running a simple script allows for automatic reconstruction, splitting, thresholding, and segmentation of stickleback specimens.
All steps in the automated pipeline are much faster than an expert operator, with minimal active time on the part of the user.
This reproducible pipeline allows for mass sampling and population-scale analysis of stickleback specimens.
Findings from Ceratobranchial Analysis
The 1st and 2nd gill rakers are remarkably different in morphology and in size and breadth, suggesting that these structures may respond differently to shifts in diet - even within the same feeding apparatus.
CB1 does have statistically significant differences as observed by the linear and pairwise analyses (and PC4, corresponding to 4.33% of the variation, appears to pull out these differences - see supplementary figures) but we caution that there are also more lakes and specimens available for this study - these differences could be due to the statistically different variances observed between the two groups.
Thankfully, this pipeline will allow for more extensive analysis from future sampling years to confirm these findings.
Within CB2, limnetic fish (and particularly those from South Rolly lake) appear to have narrower, less keeled bones than benthic fish.
The muscles that attach to the ceratobranchials (m. adductor branchialus, m. abductor filament, and m. obliquus ventralis) are attached along the side of these bones - the increased surface area in benthic fish would relate to increased muscle attachment, which would directly influence the fish to abduct these structures during water filtration (n.d.-a; Takata & Sasaki, 2005).
While dietary analyses of these fish are still ongoing, these findings suggest that the South Rolly population may have unique dietary specializations and would be expected to feed on different than fish from the other lakes in this study.
In terms of reintroduction, populations from this lake might be expected to fare better than others in terms of limnetic specialization.
Indeed, not all ecotypes represented in the FITNESS study present a uniform monolith - indeed, fish with South Rolly heritage outperform Spirit lake fish in every transplant in which they are both included (Eckert et al., 2026).
The different response of the first and second ceratobranchial brings up the possibility of a modular response to dietary shifts within the branchial basket.
Studies treating the unit as a single structure, focusing on the first ceratobranchial, or investigating external morphology might potentially miss significant changes in shape and size of these structures - and this pipeline has provided the investigators with a wealth of data with which to perform a follow-up study.
Future improvements and issues
As with many multiscan projects, the scanning parameters can be tooled individually for each scan but not for each individual specimen.
In addition, atypically large or dense specimens cause an issue for the holder and the replicability across scans.
As with most machine learning approaches, it is also important to ensure that the variation ranges across the entire dataset is represented in training to avoid erroneous segmentations.
Conclusion
The provided pipeline for the analysis of stickleback specimens provides a repeatable, high-throughput method for the analysis of 3D shape.
Although used here on exemplary stickleback specimens, the described methods can readily be applied to mass sampling efforts of multiple taxonomic groups, with data acquired on many different micro-CT machines and different sample holders, due to the use of simple region detection and reproducible logging.
The cropping out of individual specimens from a multiscan is efficient, using a custom 3D-printed sample holder and associated splitting to reduce the need for active time and maximized automated processing.
The reproducible scans and their consistent quality rapidly provide a large amount of similar data ideal for training machine learning models.
Biomedisa, as currently applied, performs on average five times faster than an skilled operator, and without the inter-operator bias endemic to splitting this amount of specimens across multiple investigators.
This brings virtual, non-destructive dissection of internal stickleback up to parity with hand-dissected methods.
Finally, the 3D analysis step of the pipeline allows for insights from 3D data that are unable to be gleaned from traditional dissection, including complex shapes and arrangements not possible under destructive sampling regimes.
We are grateful to the Microscopy Imaging Center of the University of Bern for their infrastructural support.
We also thank the manubot project (Himmelstein et al., 2019) for facilitating collaborative writing of this manuscript.
Supplementary Materials
Parameters of tomographic scans of all the fish
The CSV file ScanningDetails.csv gives a tabular overview of all the (relevant) parameters of all the scans we performed.
This file was generated with the data processing notebook and collates the relevant data read from all the log files of all the scans we performed.
A copy of each log file containing all scanning parameters is available in a folder in the data processing repository.
Adams, D. C., & Otárola‐Castillo, E. (2013). geomorph: an<scp>r</scp>package for the collection and analysis of geometric morphometric shape data. Methods in Ecology and Evolution, 4(4), 393–399. https://doi.org/10.1111/2041-210x.12035
BERNER, D., ADAMS, D. C., GRANDCHAMP, A. ‐C., & HENDRY, A. P. (2008). Natural selection drives patterns of lake–stream divergence in stickleback foraging morphology. Journal of Evolutionary Biology, 21(6), 1653–1665. https://doi.org/10.1111/j.1420-9101.2008.01583.x
Brooks, H., Haines, G. E., Lin, M. C., & Sanderson, S. L. (2018). Physical modeling of vortical cross-step flow in the American paddlefish, Polyodon spathula. PLOS ONE, 13(3), e0193874. https://doi.org/10.1371/journal.pone.0193874
Dask Development Team. (2016). Dask: Library for dynamic task scheduling. https://dask.org
Eckert, L., Bolnick, D. I., Derry, A. M., Haines, G. E., Heckley, A. M., Lind, Å. J., Peichel, C. L., Roth, A. M., Steinel, N. C., Vlahiotis, K., Weber, J. N., Hendry, A. P., & Barrett, R. D. H. (2026). Intrinsic fitness differences outweigh environmental matching in shaping introduction outcomes in nature. openRxiv. https://doi.org/10.64898/2026.02.04.699496
Ellis, N. A., & Miller, C. T. (2016). Dissection and Flat-mounting of the Threespine Stickleback Branchial Skeleton. Journal of Visualized Experiments, (111). https://doi.org/10.3791/54056
Ford, K. L., Albert, J. S., Summers, A. P., Hedrick, B. P., Schachner, E. R., Jones, A. S., Evans, K., & Chakrabarty, P. (2023). A New Era of Morphological Investigations: Reviewing Methods for Comparative Anatomical Studies. Integrative Organismal Biology, 5(1). https://doi.org/10.1093/iob/obad008
Haines, G. E., Stuart, Y. E., Hanson, D., Tasneem, T., Bolnick, D. I., Larsson, H. C. E., & Hendry, A. P. (2020). Adding the third dimension to studies of parallel evolution of morphology and function: An exploration based on parapatric lake‐stream stickleback. Ecology and Evolution, 10(23), 13297–13311. https://doi.org/10.1002/ece3.6929
Hendry, A. P., Barrett, R. D. H., Bell, A. M., Bell, M. A., Bolnick, D. I., Gotanda, K. M., Haines, G. E., Lind, Å. J., Packer, M., Peichel, C. L., Peterson, C. R., Poore, H. A., Massengill, R. L., Milligan‐McClellan, K., Steinel, N. C., Sanderson, S., Walsh, M. R., Weber, J. N., & Derry, A. M. (2024). Designing eco‐evolutionary experiments for restoration projects: Opportunities and constraints revealed during stickleback introductions. Ecology and Evolution, 14(6). https://doi.org/10.1002/ece3.11503
Himmelstein, D. S., Rubinetti, V., Slochower, D. R., Hu, D., Malladi, V. S., Greene, C. S., & Gitter, A. (2019). Open collaborative writing with Manubot. PLOS Computational Biology, 15(6), e1007128. https://doi.org/10.1371/journal.pcbi.1007128
Kikinis, R., Pieper, S. D., & Vosburgh, K. G. (2013). 3D Slicer: A Platform for Subject-Specific Image Analysis, Visualization, and Clinical Support. In Intraoperative Imaging and Image-Guided Therapy (pp. 277–289). Springer New York. https://doi.org/10.1007/978-1-4614-7657-3_19
Kluyver Thomas, Ragan-Kelley Benjamin, Pérez Fernando, Granger Brian, Bussonnier Matthias, Frederic Jonathan, Kelley Kyle, Hamrick Jessica, Grout Jason, Corlay Sylvain, Ivanov Paul, Avila Damián, Abdalla Safia, Willing Carol, & Jupyter Development Team. (2016). Jupyter Notebooks &ndash; a publishing format for reproducible computational workflows. In Positioning and Power in Academic Publishing: Players, Agents and Agendas. IOS Press. https://doi.org/10.3233/978-1-61499-649-1-87
Lösel, P. D., van de Kamp, T., Jayme, A., Ershov, A., Faragó, T., Pichler, O., Tan Jerome, N., Aadepu, N., Bremer, S., Chilingaryan, S. A., Heethoff, M., Kopmann, A., Odar, J., Schmelzle, S., Zuber, M., Wittbrodt, J., Baumbach, T., & Heuveline, V. (2020). Introducing Biomedisa as an open-source online platform for biomedical image segmentation. Nature Communications, 11(1). https://doi.org/10.1038/s41467-020-19303-w
McGee, M. D., Schluter, D., & Wainwright, P. C. (2013). Functional basis of ecological divergence in sympatric stickleback. BMC Evolutionary Biology, 13(1), 277. https://doi.org/10.1186/1471-2148-13-277
Miles, A., Kirkham, J., Durant, M., Bourbeau, J., Onalan, T., Hamman, J., Zain Patel, Shikharsg, Rocklin, M., Dussin, R., Schut, V., Andrade, E. S. D., Abernathey, R., Noyes, C., Sbalmer, Pyup.Io Bot, Tran, T., Saalfeld, S., Swaney, J., … Banihirwe, A. (2020). zarr-developers/zarr-python: v2.4.0 (Version v2.4.0) [Computer software]. Zenodo. https://doi.org/10.5281/zenodo.3773450
Overview — K3D-jupyter documentation. (n.d.). Retrieved July 1, 2026, from https://k3d-jupyter.org/
R Core Team. (2021). R: A language and environment for statistical computing. R Foundation for Statistical Computing. https://www.R-project.org/
Rawson, S. D., Maksimcuka, J., Withers, P. J., & Cartmell, S. H. (2020). X-ray computed tomography in life sciences. BMC Biology, 18(1). https://doi.org/10.1186/s12915-020-0753-2
Reid, K., Bell, M. A., & Veeramah, K. R. (2021). Threespine Stickleback: A Model System For Evolutionary Genomics. Annual Review of Genomics and Human Genetics, 22(1), 357–383. https://doi.org/10.1146/annurev-genom-111720-081402
Rolfe, S., Pieper, S., Porto, A., Diamond, K., Winchester, J., Shan, S., Kirveslahti, H., Boyer, D., Summers, A., & Maga, A. M. (2021). SlicerMorph: An open and extensible platform to retrieve, visualize and analyse 3D morphology. Methods in Ecology and Evolution, 12(10), 1816–1825. https://doi.org/10.1111/2041-210x.13669
Schindelin, J., Arganda-Carreras, I., Frise, E., Kaynig, V., Longair, M., Pietzsch, T., Preibisch, S., Rueden, C., Saalfeld, S., Schmid, B., Tinevez, J.-Y., White, D. J., Hartenstein, V., Eliceiri, K., Tomancak, P., & Cardona, A. (2012). Fiji: an open-source platform for biological-image analysis. Nature Methods, 9(7), 676–682. https://doi.org/10.1038/nmeth.2019
Schluter, D., & McPhail, J. D. (1992). Ecological Character Displacement and Speciation in Sticklebacks. The American Naturalist, 140(1), 85–108. https://doi.org/10.1086/285404
Takata, Y., & Sasaki, K. (2005). Branchial structures in the Gasterosteiformes, with special reference to myology and phylogenetic implications. Ichthyological Research, 52(1), 33–49. https://doi.org/10.1007/s10228-004-0251-5
WILLACKER, J. J., VON HIPPEL, F. A., WILTON, P. R., & WALTON, K. M. (2010). Classification of threespine stickleback along the benthic-limnetic axis. Biological Journal of the Linnean Society, 101(3), 595–608. https://doi.org/10.1111/j.1095-8312.2010.01531.x
廖炳松, 陳澤生, & 詹寶珠. (2001). A Fast Algorithm for Multilevel Thresholding. Journal of Information Science and Engineering, 17(5), 713–727. https://doi.org/10.6688/jise.2001.17.5.1
0000-0003-3388-9187
· 
habi
· 
@habi@mastodon.social
