Heritable functional architecture in human visual cortex

Highlights • We analyzed retinotopic maps from monozygotic and dizygotic twin pairs.• Visual field maps in V1-V3 are more similar in monozygotic twins.• Heritability is greater in V1 and V3 for polar angle and population receptive field sizes.• Eccentricity maps show lesser degree of heritability.• Further evidence for link between cortical morphology and topology of retinotopic maps.


Introduction
Many aspects of visual perception show pronounced individual differences.These variations have further been shown to be at least partly genetically determined in processes including binocular rivalry (Miller et al., 2010), bistable perception (Shannon et al., 2011), eye movement patterns (Constantino et al., 2017;Kennedy et al., 2017), and even complex functions like face recognition (Wilmer et al., 2010;Zhu et al., 2010).This mirrors the heritability in the coarse structural and morphological features of the brain (Chen et al., 2011).A genetic component has also been reported for aspects of cortical function thought to underlie visual processing.For instance, studies using magnetoencephalography have shown that the peak frequency of visually-induced gamma oscillations in early visual cortex is heritable (van Pelt et al., 2012).Such oscillations may derive from the local circuitry of neuronal populations (Pinotsis et al., 2013) that in turn relate to the macroscopic cortical morphology (Gregory et al., 2016;Schwarzkopf et al., 2012).Perhaps the most striking property of the human visual system is its functional organization into retinotopic maps, where adjacent locations in the visual field map onto adjacent neuronal populations in visual cortex (Glickstein and Whitteridge, 1987).This organization further determines the borders of individual regions in human visual cortex (Engel et al., 1997;Sereno et al., 1995;Wandell et al., 2007).Both intrinsic genetic and extrinsic afferent processes govern how the cortex differentiates into different regions during development (O'Leary et al., 2007;Rakic et al., 2009).Given these links, it follows that the organizational principles of visual cortex and its fine-grained functional properties may also be partially heritable.However, to date this remains untested.
Twin studies provide a unique opportunity to study environmental and genetic factors for individual differences in the architecture of the human visual system.The genetic proportion of observed variation between humans is defined as heritability (Falconer, 1965;Jansen et al., 2015), which is distinct from environmental influences, either shared or unique.Identical (MZ) twins share 100% of their genes, whereas non-identical (DZ) twins share 50% on average.Classical twin designs compare correlations between MZ and DZ twin pairs to isolate how much variance in a variable of interest can be attributed to genetic versus environmental components (Falconer, 1965).
Here, we set out to understand how much of the retinotopic organization of human visual cortex can be explained by genetic factors.We examined broad characteristics of the functional architecture of early visual cortex by conducting retinotopic mapping experiments with functional magnetic resonance imaging (fMRI) on both MZ and DZ twins.Specifically, each participant underwent three 8-minute runs of fMRI in which they viewed stimuli combining a rotating wedge with expanding or contracting rings (Alvarez et al., 2015;Stoll et al., 2020).We then used population receptive field (pRF) analysis (Dumoulin and Wandell, 2008;Moutsiana et al., 2016) to estimate the preferred visual field location and spatial selectivity of each voxel in visual cortex.We used an analysis approach with high statistical power that effectively treats each twin pair as an independent replication.We found that the organization of retinotopic maps in early regions V1-V3 is more similar in MZ than DZ twin pairs, suggesting a genetic component.

Greater similarity in retinotopic map organization for MZ than DZ twins
We first consider the heritability of broad characteristics of the retinotopic maps in early visual areas.We used spatially-normalized surface maps for each participant, generated by cortical alignment to a common template brain.Figure 1 shows example retinotopic maps from one MZ twin pair (Figure 1A) and one DZ twin pair (Figure 1B) on the occipital lobe template.The general topology of the visual regions is more similar in the MZ pair than the DZ pair, which can be seen in the overlaid maps where borders between the visual areas are more closely aligned for the MZ (Figure 1C) than the DZ (Figure 1D  To quantify the topological similarity of these visual regions, we analyzed the overlap of regions V1-V3, based on manual delineations of the retinotopic maps normalized to the template surface.To exclude experimenter bias, we carried out the delineation blinded to the genetic status of each participant.For each twin pair, we then determined the number of vertices in the reconstructed surface mesh that belonged to a given visual region in either twin.We further quantified the proportion of those vertices that overlapped between the two twins in the pair.Figure 2 shows that, on average, the overlap across all three regions was consistently greater in MZ than DZ twins (two-way analysis of variance, main effect of genotype: F(1,93)=18.29, p<0.0001).The amount of overlap also decreased from V1 to V3 (main effect of region: F(2,93)=151.22,p<0.0001).Variability in the topographical location of these visual regions is known to increase along the cortical visual hierarchy (Benson et al., 2012;Wang et al., 2015), and thus, the delineations of borders in extrastriate areas V2 and V3 are likely to be more variable than for region V1.The overlap of a particular region also depends on the borders with neighboring regions, making the results for the three regions partially dependent.Nevertheless, the difference between twin types did not differ significantly across the visual regions, with no significant interaction between the factors for visual region and twin group (F(2,93)=0.22,p=0.8035).
We further quantified the significance of these overlap statistics using a permutation analysis.We shuffled the individual participants 1,000 times and thus compared each participant to a pseudo-randomly chosen participant instead of their twin.We conducted this analysis separately for each twin group.For each visual region and twin group, we then calculated the proportion of permutations in which the overlap between pseudo-random pairs was equal to or greater than the overlap measured in the twin pairs.This effectively determines the probability of observing the level of overlap we found between twins in unrelated participants.While the overlap was significantly above this level for MZ twins in all cortical regions (V1-V3 all p<0.001), for DZ the overlap was not significantly different in any region (V1: p=0.325;V2: p=0.232;V3: p=0.429).

Heritable functional architecture in early visual regions
Our findings indicate that the general topology of visual regions is more similar in MZ than DZ twins.However, manual delineation is susceptible to errors resulting from variable signal-tonoise ratios, missing data, errors in spatial alignment, or artifacts in the pRF analysis.More importantly, this analysis only quantified the correspondence of borders between visual regions.
Next, we therefore examined the heritability of the functional architecture within each visual region, as revealed by population receptive field (pRF) analysis.For this we again used retinotopic maps spatially aligned to a common template brain.We quantified the degree of similarity in the spatial distribution of pRF parameters for all responsive vertices in the cortical surface model that surpassed an average goodness-of-fit of R 2 >0.1 in the pRF analysis across all participants.Specifically, we extracted the best-fitting polar angle and eccentricity preferences across vertices, as well as pRF sizes, separately for each visual region, in each case after removing the general trend across all participants.To this end, we first calculated an average retinotopic map across all participants irrespective of twin status, and manually delineated visual regions V1, V2 and V3 based on this grand average map.Then we used a two-stage, random-effects analysis: we first calculated the individual correlation of these pRF properties between each twin pair, for each visual region.As a measure of similarity, we calculated the circular correlation of pRF polar angle and Spearman's ρ for eccentricity and pRF size.Finally, we calculated the average correlation for each twin group and a heritability factor using Falconer's formula (Falconer, 1965).A bootstrapping test was used to determine the significance of heritability for each pRF parameter and visual region using 10,000 resamples.This analysis effectively treats each twin pair as a replication, meaning that the correlation for each twin pair is based on thousands of cortical surface vertices as observations.Therefore, we can estimate heritability with a high degree of precision.To quantify this formally, we conducted a power analysis by simulating a range of plausible true heritability levels.We simulated data for each heritability level 1,000 times, and quantified the significance of the obtained heritability estimate (Bonferroni corrected α=0.0056, for 3 pRF parameters and 3 visual regions).Then we quantified the proportion of the 1,000 simulations in which our bootstrapping test detected the heritability.Across all regions, our analysis is well conditioned with 80% power to detect a heritability of approximately 3% (Figure 3, Supplement 1).
Our results suggest a degree of heritability for all three pRF properties -their preferred polar angle and eccentricity values, and their size.Figure 3A-C shows the average intra-class correlation for the three pRF properties and the three visual regions V1-V3.For polar angle (Figure 3A) and pRF size (Figure 3C), correlations were moderate (mean r range: 0.27-0.61),suggesting a high amount of unique variance in the members of each twin pair.Correlations were nonetheless stronger for MZ than DZ twins, indicative of some heritability.In contrast, correlations for eccentricity were generally strong (Figure 3B; mean r range: 0.65-0.76)but at similar levels for MZ and DZ twins, indicative of lower heritability.
We calculated these correlations after subtracting the average pattern across all participants, which allowed us to remove general trends such as the eccentricity gradient.Nevertheless, we further used a permutation analysis to examine if any of these general patterns were driving these results.We recalculated the correlations after shuffling the participants in each group 1,000 times to break up the twin pairs.Then, we determined the significance of these correlations by quantifying the proportion of resamples in which the mean intra-class correlation in shuffled pairs was equal to or greater than the mean correlation actually observed in the twin pairs.This measures the probability of these correlations occurring in unrelated participants.Interestingly, this analysis revealed that correlations for MZ twins were significantly higher than baseline for all pRF parameters and in all regions (all p≤0.005).In contrast, the correlations for DZ twins were not significantly greater than would be expected for unrelated participants (polar angle: all p≥0.105; eccentricity: all p≥0.120; pRF size: all p≥0.435).
We further analyzed the increase in heritability across the visual hierarchy by fitting a linear regression to each bootstrap iteration of heritability values for the three regions (dummy coded as 1, 2, and 3) and determined the statistical significance of this change based on the bootstrapped slopes.This showed that the heritability for pRF size increased significantly from V1 to V3 (p=0.0047,Bonferroni corrected α=0.0167), while there was no such increase for polar angle (p=0.5283) or eccentricity (p=0.4381).
In summary, our results demonstrate that visual cortical architecture is partially heritable -both in the topographic location of visual regions, and in the broad variations of pRF parameters, particularly polar angle and pRF size.

Discussion
We examined whether the topology of early visual cortex has a genetic component by conducting fMRI testing on identical and non-identical twin pairs.Our findings suggest that both the anatomical organization and the functional architecture of early visual regions are at least partly driven by genetic components.
Our analysis revealed a moderate genetic component for the architecture of retinotopic maps, most notably for pRF sizes, a measure of the spatial selectivity of neural populations.The heritability factor for both pRF size and preferred polar angle ranged from 15-36% across regions of visual cortex.This suggests that at least two thirds of the variance in these measures was determined by environmental factors and random fluctuations.Interestingly, we observed strong correlations for pRF eccentricity preferences within twin pairs, irrespective of MZ or DZ status.This suggests only a weak heritable component of 9-21% for eccentricity gradients in retinotopic organization, suggesting that a large proportion of the variance in eccentricity preferences may be due to other factors, including experience-dependent plasticity and/or constraints set by other processes in cortical development.
In all our analyses, we used retinotopic maps that were spatially aligned to a common template brain.It is therefore likely that there are consistent patterns in retinotopic maps across all participants irrespective of twin status or familial relationship.To control for this confound, we subtracted the grand average map across the whole sample from each individual map before conducting the correlation analysis on the residuals.However, this procedure only corrected for part of the general trend.It is possible that similar visual experiences during development, such as the exposure to urban versus rural settings, could result in plasticity that shapes retinotopic map organization.However, our permutation analysis using shuffled twin pairings suggested that the distributions of all pRF parameters were much more similar between MZ twins than would be expected between unrelated participants.This was not the case for DZ twins.Therefore, we surmise that environmental factors likely only play a small role in shaping retinotopic maps.
Our findings are noteworthy in the context of previously reported links between perception and properties of retinotopic cortex.Several studies have shown that activation patterns in V1 reflect the apparent size of stimuli (specifically, the apparent stimulus eccentricity) rather than the veridical size of the retinal image (Fang et al., 2008;Murray et al., 2006;Pooresmaeili et al., 2013;Sperandio et al., 2012).Further, the macroscopic surface area of V1 correlates with variations in perceived size from phenomena like the Ebbinghaus illusion and other visuoperceptual abilities (Bergmann et al., 2015(Bergmann et al., , 2014;;Schwarzkopf et al., 2011;Schwarzkopf and Rees, 2013;Verghese et al., 2014).The rate of alternations in binocular rivalry is also heritable (Miller et al., 2010).Interestingly, V1 surface area also predicts the temporal dynamics of travelling waves in binocular rivalry (Genç et al., 2014).Given that the macroscopic surface area of visual regions correlates with the spatial selectivity (Duncan and Boynton, 2003;Harvey and Dumoulin, 2011;Song et al., 2015Song et al., , 2013)), one would expect a link between these properties and the perceived alternations in rivalry.Altogether, these findings suggest a functional role for V1 and its retinotopic organization in a range of perceptual and cognitive processes.Our finding of heritable retinotopic organization therefore suggests that there should be a heritable component for these aspects of individual perception.
The peak frequency of visually-induced gamma oscillations in visual cortex is also strongly heritable (van Pelt et al., 2012).We and others have reported links between gamma oscillation frequency, occipital levels of the neurotransmitter gamma-aminobutyric acid (Muthukumaraswamy et al., 2009), and the surface area of early visual cortex (Bergmann et al., 2016;Gregory et al., 2016;Schwarzkopf et al., 2012; but see also Cousijn et al., 2014).As a potential explanation for these links, we posited that gamma oscillation frequency depends on the cortical microarchitecture, and thus on the spatial selectivity of neuronal populations (Pinotsis et al., 2013;Schwarzkopf et al., 2012).Our present findings may constitute another missing link between these diverse aspects of visual processing: if cortical organization gives rise to gamma oscillations, then heritability in visual cortex architecture could also explain why these aspects are heritable.
The differentiation of cortex into specialized brain regions during development is driven both by intrinsic genetic mechanisms, such as signaling molecules, and extrinsic topographically organized afferents (O' Leary, 1989;O'Leary et al., 2007;Rakic, 1988;Rakic et al., 2009).In children as young as 6-7 years, the organization and functional properties of visual cortex is already similar to adults (Conner et al., 2004;Dekker et al., 2019).The development of the Stria of Gennari is preserved in congenital blindness (Trampel et al., 2011), further supporting a strong genetic determination of the arealization in V1, at least.Nevertheless, there is of course great experience-dependent plasticity in the development of functional response properties of visual cortex (Hubel and Wiesel, 1965;Wiesel andHubel, 1965a, 1965b).Early blindness causes structural changes in visual brain regions (Touj et al., 2020), possibly due to functional reorganization to serve other functions (Burton, 2003).Spatial selectivity of visual cortex, as measured by pRF size, is also abnormal in individuals with amblyopia (Clavagnier et al., 2015).Thus, it seems likely that the organization of retinotopic maps depends on both genetically determined processes and experience-dependent fine-tuning of neuronal selectivity.
Previous research has also suggested a close correspondence between cortical curvature and the structure of retinotopic maps (Rajimehr and Tootell, 2009).In fact, this correspondence allows researchers to leverage probabilistic atlases for identifying retinotopic regions without collecting functional data (Benson et al., 2012).It is possible that the constraints placed on intracortical connections due to retinotopic map organization could determine cortical folding (Van Essen, 1997).This hints at the interesting possibility that the genetic component in retinotopic map development may also drive cortical folding, though this hypothesis will need to be tested explicitly in future research.
In conclusion, our study is the first to investigate the heritability functional organization of visual regions in the human brain.Future studies must seek to understand the implications of these aspects for visual processing and perceptual function.The findings hint at an exciting possibility: how each of us perceives the visual world may be at least partly determined by genetics.

Participants
We collected data from 36 pairs of same-sex twins.All participants were healthy and had normal or corrected-to-normal visual acuity.Participants were financially compensated for their time and travel costs.We obtained written informed consent from all participants and all procedures were approved by the University College London Research Ethics Committee.
Monozygotic (MZ).We recruited 22 pairs of monozygotic (MZ) twins (18 female, 4 male), mean age 25.1 years (18-47 years).Fifteen pairs were right-handed, one pair was left-handed, six pairs were mixed (one twin right-handed, one left-handed).fMRI data from three MZ pairs was excluded: one pair was excluded due to an incidental neurological finding in one twin, with two pairs excluded due to excessive head movement during scanning.All analyses presented here used the remaining 19 MZ pairs.Dizygotic (DZ).We recruited 14 pairs of dizygotic (DZ) twins (10 female, 4 male), mean age 25.0 years (18-40 years).Thirteen pairs were right-handed, and one pair was left-handed.All DZ twins are included in the analyses presented here.

Questionnaire
Participants were also asked to complete a questionnaire to provide information about environmental factors in their upbringing, such as how often they shared school classes and friends, dressed alike as children, and how often they keep in contact with their twin.In addition to self-report of zygosity from the twins (Question 2), this questionnaire included two questions previously validated to classify twin zygosity with 95% accuracy (Sarna et al., 1978): During childhood, were you and your twin as alike as 'two peas in a pod' or were you of ordinary family likeness?(Question 3) Were you and your twin so similar in appearance at school age that people had difficulty in telling you apart?(Question 4) Demographics and questionnaire responses.Self-report of twin zygosity originally led to 17 MZ and 19 DZ pairs of twins.The questionnaire results found conflicting categorizations, where the questionnaire criteria zygosity for some twin pairs conflicted with self-reported zygosity.Five DZ pairs disagreed on both criteria (Questions 3 and 4), and three DZ and two MZ pairs disagreed on one criterion (Question 3 or 4).These 10 twin pairs were contacted and asked if they would take a genetic test through a third-party company (NorthGene, UK).Seven of the eight DZ pairs took the test, and of those, genetic testing found that the five who disagreed on both criteria were (probable) MZ twins.The two (of three) DZ pairs who disagreed on one criterion were both DZ pairs.None of the MZ twin pairs took the genetic test.We consequently switched the categories for the DZ twins whose genetic testing indicated they were MZ twins, leaving us with 22 MZ and 14 DZ twin pairs.

Parameters
Imaging data were collected on a Siemens Avanto 1.5T MRI scanner located at the Birkbeck-UCL Centre for NeuroImaging, using a 32-channel head coil with the two coils in the middle of the top half of the head coil restricting vision removed, leaving 30 channels.Functional data were acquired with a T2*-weighted multiband 2D echo-planar sequence (2.3 mm isotropic voxels, 36 slices, FOV = 96×96 voxels, TR = 1s, TE = 55ms, acceleration factor = 4).Slices were oriented to maximize coverage of occipital cortex, generally approximately parallel to the calcarine sulcus.Each participant completed three functional runs mapping population receptive fields (pRF) with 490 volumes per run (including 10 dummy volumes), and two functional runs for localizing face and scene regions (not reported here).A high-resolution T1-weighted MPRAGE structural scan (voxels = 1mm isotropic, TR = 2730ms, TE = 3.57ms) was also obtained for each participant.

Stimuli and task
Each scanning session lasted approximately 1 hour.Participants lay supine with the stimuli projected onto a screen (resolution: 1920 × 1080) at the back of the bore, via a mirror mounted on the head coil.The total viewing distance was 68cm.
We used a wedge and ring stimulus containing colorful images to map participants' visual field locations.The wedge subtended a polar angle of 12° and rotated in 60 discrete steps (one per second).The maximal eccentricity of the ring was 8.5° and expanded/contracted over 36 logarithmic steps.Within reach run, there were 6 cycles of wedge rotation and 10 cycles of ring expansion/contraction, interleaved with a 30 s fixation-only period after every quarter of the run.This stimulus rotated, expanded, and contracted around a central fixation dot.The order of rotation and expansion/contraction was the same in each run.The wedge and ring apertures contained previously described natural images (Moutsiana et al., 2016) or phase-scrambled versions thereof.Every 15 s the stimuli alternated between intact and phase-scrambled images.The sequence of individual images was pseudo-randomized.Participants were instructed to fixate at all times, and press a button if the fixation dot changed color from black to red (it could also change to a range of other colors) or if they saw a tartan pattern appear within the wedge and ring stimulus.

Pre-processing and pRF modelling
Functional MRI data were pre-processed using SPM12 (http://www.fil.ion.ucl.ac.uk/spm),MATLAB (MathWorks), and our custom SamSrf 5 toolbox for pRF mapping (https://doi.org/10.6084/m9.figshare.1344765).The first 10 volumes of each run were removed to allow the signal to reach equilibrium.Functional images were mean bias corrected, realigned and unwarped, and co-registered to the structural scan, all using default SPM12 parameters.FreeSurfer (https://surfer.nmr.mgh.harvard.edu/fswiki)was used for automatic segmentation and reconstruction to create a 3D inflated model of the cortical from the structural scan.Functional data were projected onto the reconstructed cortical surface mesh, by sampling for each mesh vertex the time course from the nearest voxel midway between the white and grey matter surface.Linear trends were removed and time courses were z-normalized.The time courses of the three pRF mapping runs were averaged.Only vertices in the occipital lobe were included for further analyses, and all further analyses were performed in surface space.
Three parameters of a symmetrical, two-dimensional Gaussian pRF model were estimated for each voxel independently: x0, yo, and σ, where the first two denote the center coordinates of the pRF in the visual field and the third is the estimate of pRF size (standard deviation of Gaussian).The model predicted the neural response at each time point of the fMRI time course from the overlap between the pRF model and a binary mask of the visual stimulus; the resulting time course was then convolved with a canonical hemodynamic response function.We then found the combination of pRF parameters whose time course best predicted the actually observed time course.Various descriptions of the data were then derived from these parameters, including: polar angle, eccentricity, and R 2 (proportion variance explained).
We conducted pRF model fitting in two stages.First, a coarse fit using an extensive grid search was performed on data smoothed with a large Gaussian kernel on the spherical surface (FWHM=5 mm).The best fitting pRF parameters were determined as those producing the maximal Pearson correlation between the predicted and observed fMRI time course.Then we conducted a fine fit, using parameters identified by the coarse fit on a vertex by vertex basis to seed an optimization algorithm (Lagarias et al., 1998;Nelder and Mead, 1965) to minimize the sum of squared residuals between the predicted and observed time course.Only vertices whose goodness of fit on the coarse fit exceeded R 2 >0.05 were included in the fine fit.This stage used the unsmoothed functional data and also included a fourth amplitude parameter to estimate response strength.The final estimated parameter maps were then again smoothed on the spherical surface (FWHM=3 mm).

Spatial normalization and manual delineation of visual regions
In order to compare retinotopic maps directly across different participants, we aligned all individual surfaces to the common space of the FreeSurfer fsaverage template.We calculated an average retinotopic map separately for each group, MZ and DZ, respectively.Then we averaged these two group maps together into one grand average map.This minimizes the undue influence the MZ group could have had on the average map due to its larger sample size.We then delineated visual regions V1, V2, and V3 based on reversals in the polar angle map and the extent of the activated portion of visual cortex along the anterior-posterior axis.Furthermore, we delineated the maps of all individual participants.This delineation was conducted blind with regard to the twin status by presenting the maps of individual MZ and DZ in a shuffled order and hiding any identifying information from the analyst.
Data in regions V3A, V3B, and V4 were less consistent across participants, and particularly susceptible to variable signal-to-noise ratios between participants.Therefore, we did not analyze data from these regions further and restricted all of our analyses to V1-V3 only.

Correlation analysis of pRF parameters
We then used the data extracted from each region to calculate the Spearman correlation of eccentricity and pRF size, and the circular correlation for polar angle, respectively, between twins in each pair.For each twin pair, we only used pRF data that surpassed an average goodnessof-fit threshold of R 2 >0.1 in the pRF model fitting across all participants.It is crucial to remove global trends before conducting such an analysis; otherwise any correlation between twins could simply be due to the pattern shared by all participants.For instance, the gradient of the eccentricity map will generally increase along the posterior-anterior axis.We therefore first removed the global trend from these data by computing the mean pattern of pRF parameters for the whole participant sample, and subtracting this mean pattern from the data of each individual.The pair-wise correlation was then carried out using the resulting residual parameters.
Next, we calculated the average correlation across each twin group (after Fisher's ztransformation to linearize the correlation coefficients) and used bootstrapping to determine the 95% confidence intervals of these averages by resampling the correlation coefficients from twin pairs 10,000 times with replacement and determining the 2.5 th and 97.5 th percentile of the resulting distribution.Heritability (H 2 ) was calculated for each bootstrap iteration using Falconer's formula (Falconer, 1965) given by H 2 = 2(M(rMZ) -M(rDZ)), where M(rMZ) and M(rDZ) stand for the mean across correlation coefficients in the MZ and DZ groups, respectively.
Finally, we conducted statistical inference on this using the bootstrapped data.To test if the heritability within a given visual region was significantly greater than zero, we quantified the proportion of bootstrap iterations where heritability was ≤0.To determine if heritability increased across the visual hierarchy, we fit a linear regression to the heritability values (regions dummy coded as 1, 2, and 3) for each bootstrap iteration.We then determined the significance of this regression by quantifying the proportion of iterations where the slope was ≤0.

Power analysis
We conducted a simulation to estimate the statistical power of the correlation analysis.To this end, we simulated data with the same number of recording sites as the smallest visual region V3 (3409 vertices) and the same sample sizes for the MZ (n=19) and DZ (n=14) groups as in our study.We varied the difference in correlation between the MZ and DZ groups to simulate a range of true heritability values.We repeated this simulation 1,000 times.For each simulated heritability value, we then quantified the number of times that our bootstrapped heritability analysis detected a significant effect (p<0.05, Bonferroni corrected by 9 comparisons).

Figure 1 .
Figure 1.A-B.Retinotopic polar angle maps for an identical MZ (A) and a non-identical DZ twin pair (B).Polar angle maps are shown on an inflated model of the cortical surface for the left occipital lobe.Cortical surfaces were normalized by aligning them to a common template.Each plot shows data from one individual.Greyscale indicates the cortical curvature, with darker patches corresponding to sulci, and lighter patches corresponding to gyri.The pseudo-color code (see insets) denotes the preferred polar angle in the visual field for a voxel at a given cortical location, as derived from population receptive field (pRF) analysis.The transparent borders show transition boundaries for visual areas V1, V2 and V3.C-D.Schematic retinotopic maps with the borders for both twins in the MZ (C) and DZ (D) twin pairs overlaid in distinct colors to allow a direct comparison.
) twins.Differences in the DZ pair are most evident in ventral occipital cortex (lower part of panel D), though there are also visible discrepancies for the MZ twins in the dorsal regions, particularly towards the peripheral visual-field representation (topright of panel C).

Figure 2 .
Figure 2. Proportion of overlap for visual regions V1-V3 between twins in each pair.Each dot denotes the results for a twin pair and cortical region.Diamonds indicate the group means.Red: MZ twins.Blue: DZ twins.

Figure 3 .
Figure 3. Intra-class circular correlations for polar angle (A), and Spearman correlations for eccentricity (B), and pRF size (C).MZ twin pair correlations are plotted against those for DZ twins.The solid line denotes the expected correlation if the variance was determined by genetic factors.The dashed line is the identity line.Error bars denote 95% confidence intervals derived through 10,000 bootstrap samples.D. Heritability for population receptive field parameters in each visual region.Data are shown for polar angle, eccentricity, and pRF size parameters, as derived from pRF analysis.Filled circles indicate the group means.The violin plot shows the bootstrap distribution for each pRF property and visual region, and the error bars denote 95% confidence intervals.Asterisks indicate significant differences at p<0.05, after Bonferroni correction for multiple comparisons.