Rapid generation of ecologically relevant behavioral novelty in experimental cichlid hybrids

Abstract The East African cichlid radiations are characterized by repeated and rapid diversification into many distinct species with different ecological specializations and by a history of hybridization events between nonsister species. Such hybridization might provide important fuel for adaptive radiation. Interspecific hybrids can have extreme trait values or novel trait combinations and such transgressive phenotypes may allow some hybrids to explore ecological niches neither of the parental species could tap into. Here, we investigate the potential of second‐generation (F2) hybrids between two generalist cichlid species from Lake Malawi to exploit a resource neither parental species is specialized on: feeding by sifting sand. Some of the F2 hybrids phenotypically resembled fish of species that are specialized on sand sifting. We combined experimental behavioral and morphometric approaches to test whether the F2 hybrids are transgressive in both morphology and behavior related to sand sifting. We then performed a quantitative trait loci (QTL) analysis using RADseq markers to investigate the genetic architecture of morphological and behavioral traits. We show that transgression is present in several morphological traits, that novel trait combinations occur, and we observe transgressive trait values in sand sifting behavior in some of the F2 hybrids. Moreover, we find QTLs for morphology and for sand sifting behavior, suggesting the existence of some loci with moderate to large effects. We demonstrate that hybridization has the potential to rapidly generate novel and ecologically relevant phenotypes that may be suited to a niche neither of the parental species occupies. Interspecific hybridization may thereby contribute to the rapid generation of ecological diversity in cichlid radiations.


| INTRODUC TI ON
Adaptive radiation describes the phenomenon of rapid diversification of a single lineage into an array of many species with ecologically varied adaptations (Losos, 2010;Schluter, 20002000). In the East African cichlid radiations, flocks of tens to hundreds of phenotypically diverse species have arisen in each of the major lakes in the region (Kocher, 2004;Seehausen, 2006). Phenotypic diversity in cichlids includes variation in behavior, morphology, color, and ecological specialization (Simakov et al., 2014;Greenwood, 1974;Salzburger, 2018;Seehausen, 1996).
The rate at which new species arise during adaptive radiation is often too high for new relevant mutations to emerge between successive speciation events (Hedrick, 2013). Hence, the high levels of heritable variation required in this process are likely to primarily stem from standing variation (Barrett & Schluter, 2008;Meier et al., 2017).
One way in which hybridization could facilitate adaptive diversification is by transgressive segregation-the occurrence of hybrid phenotypes with extreme trait values that exceed the range of trait values of both parental species combined (Rieseberg, Archer, & Wayne, 1999;Slatkin & Lande, 1994). While transgressive trait values may reduce hybrid fitness in both parental environments (Arnegard et al., 2014;Schluter, 2001), it can also allow the exploitation of a new habitat and niche that is unavailable to either parental species (Lamichhaney et al., 2017;Lexer, Welch, Raymond, & Rieseberg, 2003;Pereira, Barreto, & Burton, 2014;. If transgressive phenotypes can be genetically stabilized and if niche differentiation then leads to reproductive isolation from the parental species, transgressive novelty can result in speciation (Mallet, 2007;Rieseberg et al., 1999). Recently, theoretical work by Kagawa & Takimoto (2017) has demonstrated that in environments with many different ecological niches, hybridization by means of transgressive segregation can facilitate the process of adaptive radiation.
Alternatively, hybrid phenotypes can resemble one of the parental species, can be intermediate to both, or can combine parental traits into new trait combinations and the latter, as well as intermediate types could also allow occupation of an otherwise underutilized niche requiring intermediate or combinations of parental trait values (DeMarais, Dowling, Douglas, Minckley, & Marsh, 1992;Grant & Grant, 1996;Mallet, 2007).
Usually, transgressive segregation is due to segregation variance, sometimes due to overdominance or epistasis (Abbott et al., 2013;Rieseberg et al., 1999). Segregation variance typically results from recombination between parental species when quantitative trait loci (QTL) with antagonistic effects are present within each parental lineage and different between them, which have additive effects when recombined in the hybrid offspring (Abbott et al., 2013;Rieseberg et al., 1999;Rieseberg, Widmer, Arntz, & Burke, 2003;Seehausen, 2004). These requirements for transgressive segregation seem to be common in both plants and animals . For example, transgressive segregation has been demonstrated for skull shape in interspecific cichlid hybrids (Albertson & Kocher, 2005). Even when individual morphological traits do not segregate outside the parental range, new combinations of traits may still result in functional or mechanical transgression (Holzman & Hulsey, 2017;Parnell, Hulsey, & Streelman, 2008). Furthermore, covariation between traits can be relaxed in hybrids, which may facilitate expansion into new areas of morphospace through the novel combination of traits (Parsons, Son, & Craig Albertson, 2011;Selz, Lucek, Young, & Seehausen, 2014). While reduced covariation in morphological traits could result in the loss of adaptation in an already adapted population on the one hand, the associated release of populations from evolutionary constraints on the other hand may promote evolvability and adaptability in new and changed environments to which a population is not adapted yet (Parsons, Son, et al., 2011).
Here, we investigate the potential of second-generation (F2) hybrids between two generalist Lake Malawi cichlid species to exploit a novel ecological niche. The two parental species belong to two out of seven main lineages within the radiation (Malinsky et al., 2018), Astatotilapia calliptera and Protomelas taeniolatus. Previous studies of this cross have shown high levels of transgression in morphology in both first-(F1) and second (F2)-generation hybrids (Selz, Lucek, et al., 2014;, as well high levels of transgression in colorspace in the F1 hybrids. The latter resulted in assortative mating among F1 hybrids and parental species in mate choice experiments (Selz, Thommen, Maan, & Seehausen, 2014). Such transgression in hybrid female mating preferences favoring novel trait combinations in the hybrid males over those of the parental species was suggested to possibly facilitate the establishment of hybrid populations as incipient species, even in geographical proximity to the parental species (Selz, Thommen, et al., 2014).
The current study aimed to expand on these findings by testing whether the F2 hybrids of this cross between two generalist Malawi cichlids have the ability to exploit a new ecological niche distinct from both parents. Our specific test was motivated because some of the F2 hybrids in this cross had phenotypes that based on a qualitative impression of body shape and color, looked similar to those observed in a large clade of sand sifting species from Lake Malawi and were observed to frequently engage in sand sifting in the stock aquaria (OS personal observations). Sand sifting in search for food is an ecological specialization that has evolved multiple times independently and characterizes large clades of species in several of the large African great lake radiations of cichlid fish (Fryer & Iles, 1972;Konings, 2010).
To more thoroughly assess and characterize the occurrence of transgression and novelty in regard to sand sifting in these F2 hybrids, we took a novel approach that integrated experimental tests of transgressive performance with morphological and genetic analyses. Using a newly developed behavioral assay with two differently colored layers of sand, we assessed frequency and efficiency of sand sifting behavior in both parental species and their F2 hybrids. We used linear and geometric morphometrics to quantify variation in morphology. For trait mapping, we generated a linkage map from RAD-seq data and performed QTL mapping. To gauge our observations in the hybrids, we also tested individuals from some specialized sand sifting species in the same setup.
We demonstrate that in addition to transgression in morphology and the appearance of novel trait combinations, some of our F2 hybrids show high transgressive values in behaviors related to sand sifting, and may thus have the potential to occupy a novel feeding niche.

| The F2 hybrids and their parental species
One of the parental species, Astatotilapia calliptera, henceforth referred to as CAL, is an omnivorous habitat generalist (Konings, 2010).
Here, we used CAL from Chizumulu, an offshore island in Lake Malawi. The other parental species, Protomelas taeniolatus, henceforth referred to as TAE, is a generalist rock-dwelling species from Lake Malawi, feeding mainly on loosely attached algae and the meiofauna living among the algae (Konings, 2010). See Figure 1 for a picture of both parental species and their morphologies.
Two half-sib F1 hybrid families were obtained by crossing two CAL females with one TAE male. The six F2 hybrid families used in this study were obtained by crossing four males of one of the F1 families with five females of the same F1 family, and with one female of the other F1 family. Breeding protocols followed  and (Selz, Lucek, et al., 2014). All parental individuals used for breeding were derived from laboratory-bred stock populations. CAL individuals used for phenotyping and genotyping were also all derived from our stock population. Because only very few TAE individuals were available from our stock population, we added TAE individuals obtained from two aquarium breeders for phenotyping and genotyping.
We also included three females of specialized sand sifting species from Lake Malawi in the behavioral experiments: one Otopharynx tetrastigma female from our stock population, one Taeniolethrinops sp. female and one Lethrinops sp. female from the aquarium trade.
All fish were maintained and bred in a large recirculation aquarium system, with water temperature at 24-26°C and a 12:12 hr light/ dark cycle. The fish were fed flake food once a day and a custom mix of spinach, ground shrimps, spirulina powder, and vitamins once a week.

| Sand sifting behavior trials
We used the setup shown in Figure 2a to assess frequency of and efficiency in sand sifting in the F2 hybrids and both parental species.
The experimental procedure was slightly different for males and females. We observed spontaneous sand sifting in our stock tanks, independent of us feeding the fish. Thus, to avoid adding any confounding factors to our assay, we did not bury food in the sand for motivational purposes.
Each male fish was individually placed into one out of eight adjacent compartments (78 × 36 × 32 cm) separated by perforated F I G U R E 1 The two species used for the experimental cross. Shown are color photographs of representative male individuals of both species, Astatotilapia calliptera (CAL, top) and Protomelas taeniolatus (TAE, bottom). To the left of both: mean shape of TAE (ref) versus mean shape of CAL (target), and mean shape of CAL (ref) to mean shape of TAE (target), as plotted in geomorph (Adams & Otárola-Castillo, 2013) using the plotRefToTarget function. Shape changes are 1.5× magnified barriers to ensure water flow (see Figure 2a). Each compartment contained a glass form (35 × 18 × 6 cm), which was initially left empty and then filled with two different layers of sand after the fish were acclimatized overnight (18 ± 2 hr). The bottom layer of sand consisted of 500 ml black sand with a grain size of 0.5-0.8 mm.
The top layer consisted of 250 ml of beige sand with a grain size of ~1 mm, fully covering the black sand by a thin layer. The compartments themselves did not contain any sand and were covered from below with brown paper to prevent any light from shining through.
A picture of the undisturbed sand was taken before and after the sand sifting experiment, which lasted for 2 hr and was filmed from the top as well as from the front of each compartment. For filming from the top and for the pictures, we used GOPRO HERO cameras (V 01.09, resolution 1080-30). For filming from the front, we used a surveillance video system (Digi-Protect, ABUS group security center). After the experiment, the fish and the sand were removed. Each fish was tested three times (with five exceptions as we lost some individuals before the experiments were completed), always with a different neighbor and/or in a different compartment. Three to eight (mode = 6) fish were tested simultaneously in adjacent compartments. Before the first trial and between trials, the fish were kept in single compartments (divided by a thin mesh from compartments containing other fish) for a minimum of 2 weeks to give the fish time to recover between trials and to minimize potential learning effects.
Both the test and the single compartments were located in tanks that are part of the same large recirculation aquarium system where the fish were maintained (see above). Hence, the fish always had visual and olfactory contact with each other and/or with other fish.
A total of 25 F2 hybrid males, 8 CAL, and 9 TAE males were tested with this setup.
The females showed signs of discomfort when tested individually. Hence, we adjusted the protocol as follows: Instead of just one, we placed three fish of one genotype class (CAL, TAE, F2 hybrids; henceforth "class") together in a compartment, which allowed them to shoal (Pitcher, 1982), as we had observed them usually do in the stock tanks. To be able to distinguish the individuals visually, we sedated the females with MS222 (~75 mg/L, plus ~0.5 g/L sodium bicarbonate as buffer) and marked them with the ablation of a small part of fin (in the soft ray part either on the dorsal fin, on the upper or on the lower caudal fin lappet) one day before the first trial. No behavioral differences were observed between the marked and untreated females from the same tank. The compartments already contained the sand layers when we added the fish. The sand remained covered with a plastic mesh for an acclimatization time of three hours. Then we lifted the cover, took a picture of the sand, and started filming for 90 min from the front using the GOPRO HERO cameras. Afterward, another picture of the sand was taken and the fish were removed. Each group was tested twice, with a break of one day in between. 16 F2 hybrid groups, six CAL groups, and five TAE groups were tested with this setup. Two to six such groups were tested simultaneously in adjacent compartments. As activity visibly decreased from first to second trial, only the first trial per group was included in further analyses.
For a qualitative comparison with specialized sand sifting species from Lake Malawi, we also tested one group composed of one Otopharynx tetrastigma female, one Taeniolethrinops sp. female, and one Lethrinops sp. female in the same setup.

| Analysis of male sand sifting trials
For each male trial, we measured and scored six different behaviors (Table 1). All behaviors were scored by the same person (AFF). Trials with "no strike" counts combined with a "hiding" or an "interaction" score of 5 (see Table 1) were excluded from analyses (10 out of 118 trials).
The pictures of the sand tray before and after each trial were imported into ImageJ I.49v (Rasband, 2015) as a sequence and then cropped so the pictures only contained the sand area. We color-thresholded the pictures using the default thresholding F I G U R E 2 Setup for sand sifting behavior trials. (a) A tank that was part of a large recirculation aquarium system was divided up into eight adjacent compartments (only two shown) by metal plates with approx. 1 × 1 cm holes to ensure water flow. The compartments (78 × 36 cm) had a water depth of 32 cm and each contained a glass form (35 × 18 × 6 cm) which was filled with two layers of sand of different color. Where sand sifting occurred, black sand from the black bottom layer appeared on the beige sand surface (indicated as black spots on the right). Efficiency in sand sifting was calculated as the area of sand that had turned black divided by the number of strikes (how many times a fish picked up sand). (b) Visualization of color thresholding of sand pictures in ImageJ to measure the % of sand area that had been turned over-and hence turned black-in sand sifting trials. The same sand surface is shown before a trial (left) and after a trial without and with color threshold (right) method. The brightness level in the color threshold option was set to 80, which was the level we determined as best capturing the contrast between the black and beige sand particles (see Figure 2b). We then used the "analyse particles function" to calculate the percentage of area covered by black sand, that is, the area where black sand had been brought to the surface by the fish picking at or sifting through the sand (see Figure 2b). The score from the before-trial picture was subtracted from the after-trial picture. To remove effects of "scraping" (see Table 1), the score was multiplied by the corresponding inverse value of the effect due to scraping (see Table 1).
As a measure of efficiency, we divided the percentage of sand area that had been turned over by the number of strikes per individual (+0.001, to allow log-transformation). This could only be calculated for trials in which at least one strike had occurred (73 out of 108 trials).
In the classical view, transgressive values in the F2 hybrids are defined as values outside the combined parental range (Parsons, Cooper, & Albertson, 2011;Rieseberg et al., 1999;). In addition to this classical definition of transgression, we use more conservative transgression thresholds following (Parnell, Hulsey, & Streelman, 2012). Thus, we calculated transgression thresholds by taking the higher mean parental value plus two standard deviations (SD) for the high threshold, and the lower mean parental value minus two SDs for the low threshold. This approach also accounts for unequal sample sizes, as the lower sample size of the parental species should increase the SD and with it the thresholds for transgression (in both directions).
To do pairwise comparisons of means in strike counts and efficiency between the three classes (note, however, that population means are not predicted to change when transgressive variation occurs), we used the permTS function implemented in the perm package in R, with two-sided tests using the exact.mc method and a conf.
For the measure of efficiency, this analysis was done for the raw (log10) values as well as for size-corrected values given by the residuals of a linear regression (lm function) of (log10) efficiency versus (log10) standard length (SL) within each group to account for effects of size.

| Analysis of female sand sifting trials
For every individual, we recorded latency and we counted the number of strikes and scraping behavior (as defined in Table 1). Black sand area and efficiency were calculated as described above, but using the summed-up number of strikes of all three individuals in a group to calculate efficiency. Hence, the efficiency measure for females applies to the group and not to a single individual. This was possible for all groups except one group of TAE, which had zero strikes.
Additionally, the level of aggression within each group was qualitatively assessed (scores 1-5). Calculation of transgression thresholds and statistical comparisons of means were performed as outlined above.

| Linear and geometric morphometrics
All fish used in experimental trials (n = 123) were anesthetized and subsequently euthanized in MS222 (~75 mg/L for anesthesia, ~200 mg/L for euthanizing, and ~0.5 g/L sodium bicarbonate added as buffer) and subsequently stored in 75% ethanol. Another 88 F2 hybrid individuals that had died of natural causes in the aquarium between the years 2010-2016 and that had been preserved in 75% ethanol were also available for morphological analysis. A standardized photograph was taken from the left side of each fish (161 F2 hybrids, 26 CAL, 24 TAE in total; including all individuals used in the behavioral trials) using a CANON E05 60D camera with a CANON 50 mm macro lens.  Figure 3). The landmarks were set in such a way that they would be largely insensitive to slight head elevation due to lower jaw rotation, which we frequently encountered in the fish that had died of natural causes. The landmarks were used for geometric morphometric (shape) analysis and for calculating four linear distances as well as mouth angle (MA; Figure 3). The linear distances were Epaxial Depth (EpD), Eye Length (EyL), Lower Jaw Length (LJL), and Maxilla Length (ML) (Figure 3). Standard Length (SL) and

Scored behavior Description
Head Width (HW) were directly measured on the fish (average of three measurements with an accepted error of max. 5%). These traits were chosen to capture aspects of functional morphology associated with the generation of high suction pressure (McGee, Schluter, & Wainwright, 2013) without the need to clear and stain the specimens, which was not possible due to constraints on specimen use for other projects. All analysis of morphological data was performed in R, using the R package geomorph (Adams & Otárola-Castillo, 2013 To remove an effect of open mouth (see above) in geometric morphometric analysis, we rotated landmark 4 such that the angle between ML and LJL was the same for all specimen (mean angle over all specimens; (Adams, 1999). To obtain shape data, the landmark configurations were subjected to a generalized Procrustes analysis (GPA). During superimposition, the semi-landmarks (see Figure 3) were slid to minimize bending energy (Bookstein, 1997). As shape variation was significantly associated with SL due to allometric effects (Procrustes ANOVA: F 1,209 = 41.962, p = .001**), we took the residuals of a linear regression of shape coordinates on (log 10 ) SL, and added them to the consensus coordinates to obtain allometry-free shapes (Adams & Otárola-Castillo, 2013;Klingenberg, 2016).
Homogeneity of slopes of the F2 hybrids and both parental species was given (ANOVA with randomized residual permutation procedure: F 205 = 0.716, p = .77), hence size correction was performed over all three classes (F2 hybrids and both parental species) simultaneously. The major axes of shape variation were then assessed using principal component analyses (PCA) and plotted as morphospace with warpgrids (visualizing the difference between mean shape and shape at the maximum and minimum end, respectively, of a PC axis) for visual inspection of shape change along PC axes. We calculated transgression thresholds for each axis as described above (following (Parnell et al., 2012)).
Additionally, we performed the PCA only including the F2 hybrids-these PC scores were subsequently used to test for associations with efficiency in males and in QTL mapping-and we then used the predict function to project the parental species into this morphospace.
For shape comparisons with specialized sand sifting species, we placed the same set of landmarks on photographs of 21 preserved specimens (of samples stored at the Natural History Museum (atan2((Y6-Y5),(X6-X5))*(180/π))-(atan2((Y3-Y10),(X3-X10))*(180/π)). An upturned mouth will have a larger angle than a downturned mouth correction and Procrustes analysis (GPA) were performed on the whole extended dataset, size correction for these individuals was performed separately as they were generally larger and had significantly different allometry slopes (F 224 = 0.02, p = .01*). We then used the predict function to project the sand sifting species into the morphospace occupied by the F2 hybrids and both parental species.

| Associations between efficiency and morphology in F2 hybrid males
Associations of sand sifting efficiency with morphological traits in the F2 hybrid males were tested using linear models (lm function in

| RAD sequencing
Finclips of all individuals used in behavioral and morphological analyses (see overview on Datadryad) for RAD-tag sequencing had been stored in 98% ethanol. DNA was extracted from finclips using a phenol-chloroform protocol (Sambrook & Russell, 2001 Phred quality score below 10, and reads with more than 5% of bases with a quality score of less than 30. The reads of each individual were aligned to an anchored version of the Pundamilia nyererei reference genome (Feulner, Schwarzer, Haesler, Meier, & Seehausen, 2018).
End-to-end alignment was performed with Bowtie2 (Langmead & Salzberg, 2012) using default parameters. Base score recalibration was performed as described in (Marques et al., 2016

| Linkage map and QTL mapping
As finclips of the grandparents of our cross were not available, we QTL mapping was performed using the scanone function (standard interval mapping using the method EM algorithm (Lander & Botstein, 1989) in the R/qtl package (Broman, Wu, Sen, & Churchill, 2003)). Conditional genotype probabilities were calculated using the calc.genoprob function with a step-size of 1 cM, an error.prob of 0.05 and the Kosambi mapping function. Mapping was performed on all five linear traits and MA together, and on shape PC axes 1-6 together (axes explaining min. 5% of variance) using normal models. Sand sifting behaviors (maximum observed efficiency and number of strikes in males and number of strikes in females) were mapped separately using both normal and nonparametric models. Genome-wide significance thresholds were determined by permutation (n = 1,000). This was done for all fish together without and with sex as covariate for the morphological traits. Bayesian credible intervals were calculated using the baysint function. Percent variance (PVE) was calculated as 1-10 −(2/n)*LOD following (Broman et al., 2003), where n is the number of individuals.

| PERMITS
Fish experimentation and euthanasia were authorized by the veterinary office of the canton of Lucerne (License number: LU04/07).

| Transgression in sand sifting behaviors in some F2 hybrids
When analyzing the number of strikes (i.e., how often sand was picked up), we found that some F2 hybrid males had transgressive values (in the positive direction). This outcome was similar when considering the maximum observed number of strikes per individual (Figure 4a), or when conservatively considering the mean number of strikes over all three trials, or the number of strikes per trial separately ( Figure S1).
One or two males of the CAL parental species also had high values, but they were not as high as the most extreme F2 hybrids. Permutation tests to compare the means of the F2 males and the males of the two parental species revealed no significant difference between any of the three classes in the maximum observed number of strikes (F2 vs. CAL p = .754, F2 vs. TAE p = .426, TAE vs. CAL p = .104).

F I G U R E 4
Some F2 hybrids have transgressive behavioral scores. (a, b) show strike counts (i.e., how often sand was picked up by an individual) in sand sifting trials. For males (a), the maximum observed number of strikes for each individual is shown (for per trial strikes and means over all three trials see Figure S1). For females, (b), the number of strikes for each individual in the first trial is shown. (c, d) show sand sifting efficiency (log10-transformed) plotted against standard length (SL; log10-transformed). For males (c), the maximum observed efficiency for each individual is shown (per trial efficiencies and means over all three trials shown in Figure S1). For females (d), efficiency per group in the first trial is shown. (SL here is also the mean of the three females in a group). Additionally, one group including three females of three different specialized sand sifting species is shown ("Sandsifters"). The fish shown on the right are male representatives of each of the three classes, and an Otopharynx tetrastigma male is shown as representative of one of the tested specialized sand sifting species. The black dotted lines in all plots indicate transgression thresholds given by the highest/lowest parental species mean +/-2 standard deviations.  The mixed group of sand sifting species had a higher overall strike count (dashed brown line in Figure 4b) and as expected, they had a higher sand sifting efficiency score than any of the other female groups (Figure 4d).

| Transgression in several linear traits
Transgression in the F2 hybrids was observed in Epaxial Depth (EpD; both directions), Eye Length (EyL; positive direction), Mouth Angle (MA; both directions), and Head Width (HW; both directions) ( Figure 5, Table S1). However, for HW, this only concerned one individual each, and for this trait, the transgression thresholds were still within the combined parental species' range. There were some differences in the directionality of transgression when calculating the thresholds for males and females separately (see Table S1). Maxilla Length (ML) and to some extent Lower Jaw Length (LJL) differentiate the two parental species ( Figure 5): CAL have larger values in both traits.
In the PCA on all linear traits and MA ( Figure S2), some individuals fell outside transgression thresholds on several of the PC axes.
Furthermore, the F2 hybrid multivariate morphospace (as 95% confidence interval ellipses) covered a larger area than either parental species', and also covered areas outside both parental morphospaces on PC1 versus PC2, for example lower left and upper right areas ( Figure S2a; and see Table S2 for trait loadings).

| Areas of novelty in shape
In the analyses including the F2 hybrids and both parental species, PC1 and to some extent PC2 differentiated the parental species ( Figure 6). Shape change along PC1 was mainly driven by mouth size and epaxial depth, such that TAE had a smaller mouth but a deeper epaxial muscle than CAL. Shape change along PC2 was mainly driven by eye size and the angle of the mouth, such that TAE had larger eyes and a more downturned mouth. As in the PCA of linear traits, some F2 hybrid individuals fell outside transgression thresholds on several PC axes ( Figure 6). Also, the F2 hybrid morphospace (as 95% confidence interval ellipse) covered an overall larger area as well as areas outside both parental species morphospace on PC1 versus PC2 (Figure 6a). The most extreme transgression seemed to occur on PC4 (Figure 6c), where shape variation was mainly driven by mouth angle and mouth size: Individuals on the maximum end of this axis had a larger downturned mouth, while individuals on the minimum end of this axis had a smaller upturned mouth (Figure 6c). These are novel combinations compared to the parental species as CAL have a larger and more upturned mouth, and smaller eyes and a smaller epaxial depth compared to TAE (Figure 1).
When the morphospace of the 21 sand sifting species we measured was projected into the F2 hybrid and parental species morphospace, it encompassed the parental species morphospace as well as the transgressive morphospace of the F2s on PC1 versus PC2 (Figure 6b). This was mostly true also for PC3 ( Figure 6d) and PC5 versus PC6 (Figure 6f), but on PC4 (Figure 6d), a number of F2 hybrids came to lie well outside the sand sifting species morphospace.
When shape morphospace was defined by the variance in the F2 hybrids and the individuals of the parental species were predicted into this morphospace ( Figure S3; note that these were the
PC values later used in QTL mapping), there was again a substantial area where the F2 hybrids occupied morphospace not covered by either parental species (lower right corner in Figure S3a).
Compared to the minimum end (i.e., negative sign) of PC1, shape on the maximum end (i.e., positive sign) of PC1 was characterized by a large mouth (CAL-like), a large eye (TAE-like), small epaxial depth (CAL-like), and a convex dorsal head profile (more TAE-like).
On PC2, shape on the minimum end was characterized by a rather shallow head with a long snout and a small downturned mouth (overall more TAE-like). Hence, F2 hybrids in this morphospace outside both parental species morphospaces combined traits in a novel way.

| No significant associations of morphology with efficiency in F2 hybrid males
None of the linear traits or shapes were significantly associated with maximum scored sand sifting efficiency in males (all p > .1 after FDR correction; see Table S3). The best (but not significant) models with multiple traits based on AIC were a model containing shape PC axes 4-6, in which PC4 and PC6 had only marginally significant effects (p = .05 and p = .06 before FDR correction), and a model containing EyL, LJL, and ML, in which ML (a nontransgressive trait) had an only marginally significant effect (p = .08 before FDR correction) (see Table S3). Also, the behaviorally tested F2 hybrid males did not cluster by efficiency category (low, medium, high mean efficiency) in the shape PCA shown in Figure 6b,d,f.

| A marginally significant QTL for sand sifting
We found one marginally significant QTL for the number of strikes in males on LG14 (Pun-LG22/Ore-LG12 (p = .072; Table 2, Figure 7) using the nonparametric model in mapping. Homozygotes for the CAL allele had the highest scores, homozygotes for TAE the lowest.

| Significant QTLs for Epaxial Depth (EpD) and shape PC2
We found one significant QTL for EpD (p = .021; Table 2, Figure 8) explaining 29.4% of variance on LG17 (Pun-LG6, Ore-LG11) when accounting for sex as additive and interactive covariate (allows average phenotype as well as QTL effect to be different between the two sexes). At the most significant marker for this QTL, homozygotes for the CAL allele had a smaller EpD than homozygotes for the TAE, and heterozygotes had the largest values (Figure 8b).
We found two significant QTLs for shape PC2: one on LG22 (Pun-LG14, Ore-LG9) explaining 16.6% of variance in the mapping without sex as covariate, and one on LG13 (Pun-LG20, Ore-LG4) explaining 20.1% of variance, when accounting for sex as additive and interactive covariate (p = .038; Table 2, Figure 8). At the most significant markers for these QTLs, homozygotes for the CAL allele had a low PC2 score and homozygotes for the TAE allele had a high PC2 score in both cases (Figure 8). We found another marginally significant (p < .1) QTL for PC2 on LG10 (Pun-LG17, Ore-LG20) explaining 14.6% of variance, either when mapping without covariates, (p = .067), or when accounting for sex as additive covariate (p = .077) (Table 2, Figure 8). At the most significant markers for this QTL, homozygotes for the CAL allele had a high PC2 score and homozygotes for the TAE allele had a low PC2 score.

| D ISCUSS I ON
We experimentally investigated transgressive segregation in morphology and behavior that may result in the potential to exploit a novel ecological niche (i.e., sand sifting) in a second-generation (F2) hybrid cross between two trophic generalist species of Lake Malawi cichlids, in which some of the hybrids phenotypically (in body shape and color) resembled sand sifting species in the adaptive radiation of Lake Malawi cichlids (Fryer & Iles, 1972;Konings, 2010 Konings, & Sato, 2001). Interestingly, most of the specialized sand sifters among Lake Malawi cichlids are also male bower builders (Fryer & Iles, 1972;York, Patil, Hulsey, Streelman, & Fernald, 2015).
So evolutionary interactions between the two functions (feeding and reproduction) may have played a role in the adaptive radiation.
Bower building has not been observed in either parental species of our cross. We think we can exclude a (transgressive) motivation for bower building to explain our data because such behavior would only be expected in reproductively active males and usually only well after they have established a territory (Fryer & Iles, 1972;York et al., 2015).
Second, we found individuals among the hybrids with transgressive values in several linear morphological traits. As expected from theory ) and demonstrated in cichlids  combinations, which may provide the raw material for rapid adaptive diversification (Grant & Grant, 1992;Lewontin & Birch, 1966;Stebbins, 1959).
The absence of any simple significant association between morphology and sand sifting efficiency likely reflects the complexity of the functional relationships associated with the behavior (Kane & Higham, 2015). The functional basis of benthic foraging on single targeted prey items is relatively well understood (McGee et al., 2013), but the factors governing more indiscriminate sifting of large quantities of sand are more complex. Suction potential can be enhanced via enlarged epaxial musculature or a reduction in buccal volume (McGee et al., 2013), which can be seen in one of our parental species, P. taeniolatus (TAE), which has both a small gape and enlarged epaxial musculature relative to A. calliptera (CAL) (see Figure 1). However, while a smaller buccal volume may enhance suction, it also limits the amount of substrate processed per strike, compromising efficiency. This suggests that sand sifters may experience a trade-off between suction potential and sifting efficiency, producing an adaptive ridge in morphospace where no one trait combination is optimal. Furthermore, other behavioral aspects that we did not measure in our experiment, such as at which angle a fish approaches the sand, may influence efficiency in sand sifting.
We found significant QTLs for one linear trait (Epaxial Depth, EpD) and for shape PC2. The QTL for EpD could to some extent be predicted to generate transgressive phenotypes via overdominance since heterozygous genotypes had the highest values in this trait (Figure 8b). However, although three individuals were indeed transgressive in this direction for this trait (Figure 5a), more individuals were transgressive in the opposite direction, and these low values occurred in homozygous CAL genotypes, indicating that other mechanisms, such as epistasis, could be involved in generating transgression in this trait. On shape PC2, the more deep-headed individuals with a large upturned mouth and short snout (more CAL-like) were on the maximum end of PC2, the more slender-headed individuals with a small downturned mouth and a longer snout (more TAE-like) were on the minimum end of PC2 (see Figure S3a). The genotypes at the marginally significant QTL on LG10 (Pun-LG17/Ore-LG20) are consistent with this: individuals that are homozygous for CAL alleles (AA in Figure 8b) Figure 4 and Figure   S1).
In conclusion, we demonstrate a case where hybridization between two Lake Malawi cichlid species that are trophically unspecialized has not only resulted in transgression in morphology, but also in transgressive trait values in behaviors related to a feeding behavior that is not typical for either parental species but characterizes a whole clade of trophically specialized species in the Lake Malawi cichlid radiation. We find several QTLs for morphology, suggesting the existence of some loci with moderate to large effects. The genetic basis of transgression in this cross, however, is not yet resolved. One issue in our study was the relatively small sample sizes, which may make detection of QTLs difficult, especially those with small effect (see e.g., Slate, 2017), and will usually lead to overestimation of the effect of significant QTLs (Beavis, 1995). Furthermore, in addition F I G U R E 8 QTL mapping suggests the presence of some loci with moderate to large effects on morphology. In (a) LOD scores across the 22 chromosomes are shown for EpD and PC2 (mapping with sex as additive (add)/interactive(int) covariates for EpD and PC2, and without covariates (no cov) for PC2), for which we found significant (p < .05) or marginally significant (p < .1) QTLs. The dotted lines represent genome-wide significance thresholds of p < .05. (b) shows the phenotypic distributions and effect plots of these two traits at the most significant marker associated with a QTL. Genotype AA corresponds to homozygous for CAL alleles, BB to homozygous for TAE alleles, AB to heterozygous LG17/Pun-LG6/Ore-LG11 LG13/Pun-LG20/Ore-LG4 L G22/Pun-LG14/Ore-LG9 LG10/Pun-LG17/Ore-LG20 to complementary gene action and overdominance, epistasis could also be involved in generating transgression (Rieseberg et al., 1999).
The latter is still challenging to infer and would need more elaborate analyses (see e.g., Laurie, Wang, Carlini-Garcia, & Zeng, 2014).
Nevertheless, the functional morphology and/or behavior, which may underlie the effective and efficient use of a resource neither parental species has specialized on, is present in some of our F2 hybrids. If this occurred in the wild, in the situation where specialized sand sifters were not yet present, it may permit some hybrid genotypes to tap into niches not used by either parental taxon, as has also been demonstrated for natural hybrids in several plant taxa (Anton, R. Ward, & Cruzan, 2013;Donovan, Rosenthal, Sanchez-Velenosi, Rieseberg, & Ludwig, 2010;Xing et al., 2014;Zhao et al., 2014) and for experimental hybrids in Lake Victoria cichlid fish . If hybrids have elevated fitness in some ecological conditions, but are less fit than the parentals in the latter's ecological conditions, this might enable ecological niche partitioning between a population of some hybrid genotypes and both parental species (Seehausen, 2004). If niche differentiation leads to habitat isolation, this could be one route to the evolution of reproductive isolation between a hybrid lineage and its progenitors (Buerkle, Morris, Asmussen, & Rieseberg, 2000;Mallet, 2007;Rieseberg et al., 1999).
The alternative route to hybrid speciation would be if niche differentiation becomes genetically coupled to behavioral reproductive isolation (Selz, Thommen, et al., 2014). Interestingly, we see evidence for both in the hybrid cross we studied here: we see transgression and novelty in these hybrids in morphology (Selz, Lucek, et al., 2014;; and the present study) and in a behavior linked to utilizing a resource not typical of either parent (present study), and also in traits that are relevant for reproductive behavior (assortative mate choice and transgressive male color in F1 hybrids in Selz, Thommen, et al., 2014). If both assortative mating and the potential to adapt to a novel food resource come together in hybrids-and indications for both exist in this experimental species cross-this might enable a new hybrid species to emerge relatively quickly. Together, these findings suggest that interspecific hybridization in cichlids can generate functional novelty that is known from other cichlid species to be adaptive in certain ecological contexts and could thereby possibly facilitate or have facilitated niche shifts, promoting speciation, and adaptive radiation.

ACK N OWLED G EM ENTS
Sequencing was done at the Center of Integrative Genomics (CIG), University of Lausanne. RAD-seq data were processed and analyzed on clusters managed by the Genetic Diversity Centre (GDC) at ETH Zurich. Many thanks to the FishEc group for many valuable discussions and help throughout this project. We further thank Katie Peichel and two anonymous reviewers for their feedback on previous versions of this manuscript. This project was funded by SNF grant 1535 17-870 to OS.

CO N FLI C T O F I NTE R E S T
We declare we have no competing interests.