Introduction

The need for neural systems to be robust on one hand but have plasticity to adapt to experience on the other, places opposing pressures on the formation of the brain’s neural architecture. Historically, there have been two central contrasting views: early studies suggested that the brain exhibits remarkable plasticity, with the capacity for significant rewiring and reorganisation in response to new or altered experiences (Baker et al., 2005, 2008; Baseler et al., 2002; Calford et al., 2005; Dilks et al., 2009; Ferreira et al., 2016; Schumacher et al., 2008). More recent work portrays the brain as a fundamentally stable structure resistant to functional reorganization even when primary sensory inputs are lost (Adams et al., 2007; Baseler et al., 2011; Goesaert et al., 2014; Molz et al., 2023; Ritter et al., 2019). Our study reframes this debate by investigating how adaptations of higher-order processing stages can compensate for stability of sensor-to-cortex connections, offering a new perspective on neuroplasticity through the lens of hierarchical network dynamics.

In the visual domain, the focal point of the debate on plasticity and stability has hinged on the extent to which retinal input deprivation can drive local reorganisation in early visual cortex (Adams et al., 2007; Baker et al., 2005, 2008; Baseler et al., 2002, 2011; Calford et al., 2005; Dilks et al., 2009; Dumoulin & Knapen, 2018; Ferreira et al., 2016; Goesaert et al., 2014; Haak et al., 2015; Molz et al., 2023; Ritter et al., 2019; Schumacher et al., 2008) . In reality visual impairment is a more global phenomenon with complex dynamics beyond constricted local retinocortical projection zones (Carvalho et al., 2019). Investigating these dynamics is crucial for understanding neural plasticity bottlenecks of regenerative treatments under trial in human patients (Fischer et al., 2020; Georgiou et al., 2021; Michaelides et al., 2023). In these treatments, signal is recovered at the level of the sensory organ, here the eye. If tissue at any downstream processing stage has deteriorated due to inactivation or has adopted new roles, it may be less responsive to its original signals. For ground-breaking sight recovery treatments to reach their full potential, understanding these cascading dynamics in the human brain is now paramount.

We revisit the debate on neuroplasticity and stability by investigating adaptations to altered sensory input across hierarchical cortical processing stages, using the inherited retinal disease achromatopsia as a model. Achromatopsia causes cone photoreceptors in the retina to be inactive from birth (Aboshiha et al., 2014). This means visual development is purely driven by rod photoreceptors. Rod-mediated vision causes low acuity, a blind spot (scotoma) in the rod-free retina of the fovea, and no colour vision. In most studies on how the brain adapts to sight loss, a major challenge is accurately matching the sensory input experienced by patients in control individuals with typical sight (Brown et al., 2021; Jones et al., 2020;

Macnamara et al., 2021). This matching is critical for distinguishing neural adaptations to long-term altered visual experience from normal transient responses to varying visual input (Binda et al., 2013; Carvalho et al., 2021; Haak et al., 2012; Hummer et al., 2018). While normally sighted individuals rely on cone vision most of their waking time, rod only vision can also naturally occur or be induced experimentally (e.g. Barton & Brewer, 2015). This makes it possible to directly compare lifelong rod-only vision in achromatopsia with transient rod- only experiences in normal sight, distinguishing the long-term neural impacts of altered sensory input from temporary changes within the typical system’s range of flexibility.

It remains unclear how the brain adapts to altered retinal input in achromatopsia. Initially, a functional MRI study showed that areas of the early visual cortex, typically receiving only foveal cone input, might reorganise to process input from peripheral rod-containing regions (Baseler et al., 2002). This suggests potential cortical plasticity. However, recent findings by Molz et al. (2023) present a conflicting view, showing no reorganization in the deprived primary visual cortex of achromats, and structural signs of deprivation (Molz et al., 2022).

This suggests functional stability. Here, we reasoned that both stability and plasticity in a deprived visual system can co-exist across different cortical processing stages; A high degree of stability in sensor-cortex connections (Haak et al., 2015; Molz et al., 2023) may be offset by plasticity in downstream regions in the hierarchy, which may exhibit more scope to flexibly adapt to optimise the use of available input.

To test this prediction, we drew on a uniquely large sample of achromats and controls. We first characterised spatial tuning and remapping of V1 in a cone-input deprived visual system to resolve ambiguity of earlier findings. We then used a connective field modelling approach to assess sampling of V1 information by higher order visual cortex (V3). To truly distinguish long-term adaptations from transient changes to altered input, we took a different approach than prior studies. Rather than only manipulating luminance levels, we utilised the method of silent-substitution (Estévez & Spekreijse, 1982) to create photoreceptor-specific population receptive field (pRF) mapping that only stimulates rods. Crucially, this allowed us to obtain well-matched, high-quality rod-driven maps from both achromats and controls, enabling unbiased comparison of rod-mediated responses across the groups.

Using these specialised stimuli to characterise the impacts of purely rod-based vision on the visual systems of achromats and normally-sighted controls, we first show that V1 retinotopy presents the hallmarks of a highly robust system. There is no evident remapping of cortical regions that are input-deprived under rod-only viewing circumstances, despite their life-long altered input in achromatopsia. We found that the robustness to long-term altered input in V1, is paired with reorganisation of retinotopic read-out by higher-order visual cortex (V3). A comparison of the visual space organization of V1’s connectivity with V3, reveals a shift in the distribution of resources in achromats: more cortical resources in V3 are dedicated to sampling the V1 cortical surface that encodes 1°-4° visual angles - the most foveal areas accessible to rods. V3 also samples from fewer V1 vertices, effectively increasing its spatial resolution. These changes in the rod-based visual system of achromats could act to compensate for the absence of two features of the cone system that the rod system lacks: its high-density foveal sampling and its overall higher resolution. This may reflect an overarching broader organising principle where hierarchical optimisation processes allow the brain to make more efficient use of the available (in this case atypical) input.

Results

In three results sections, we characterise the interplay between neural plasticity and stability in achromatopsia across three complementary neural markers: brain structure, spatial tuning characteristics of the functional MRI signal in V1, and the hierarchical retinotopic connectivity patterns that indicate how higher-order visual cortex (V3) sample information from V1.

Structural analysis

As a proxy of structural integrity, we compared cortical thickness across 8 eccentricity segments of V1 (0°-1°,1°-2°,2°-3°, etc.; see Figure 1B) in a large sample of 31 achromats and 100 controls aged 6-35 years. Data were collated over two scanners. We expected changes in achromats to be largest in the central 0°-2° eccentricities of the visual field, encompassing the retinal ‘rod-free’ zone, an oval of approximately 1.25°x0.96° in size (Curcio et al., 1990). V1 segments were defined using the Benson atlas (Benson & Winawer, 2018), which is fitted using cortical landmarks on individual brains. The atlas fit produces predicted eccentricities to each vertex, this allows us to identify where the representation of these eccentricities would have been given typical development and investigate what they represent in atypical development.

Figure 1

Figure 1A shows an increase in cortical thickness in achromats compared to controls only in the foveal projection zone (Controlling for scanner, overall cortical thickness, and age; 0°-1° ANCOVA group effect: F_(1,126)=9.05, p_unc=0.003; 1°-2°: F_(1,126)=7.79, p_unc=0.006). This increase was not present in more eccentric parts of V1 (p_unc>0.2 for all). This pattern is in line with previously reported V1 foveal thickening in adult achromats (Molz et al., 2022), and moreover suggests that this change is already present in late childhood. A further comparison of the cortical thickness in the foveal projection zones (0°-2°) in V2 and V3 (Figure 1A right), revealed a trend toward significance in V2 (F_(1,26)=3.07, p=0.08) and a significant effect in V3 (F_(1,126)=6.32, p=0.01). As thicker cortex in early visual areas is a hallmark of visual deprivation from birth (Jiang et al., 2009; Park et al., 2009), this indicates that the foveal projection zones in achromats may lack retinal input. This provides evidence for stability, because if these cortical areas had substantially reorganised to receive more peripheral retinal input, we would not expect to see this marker of deprivation.

pRF properties of V1

We next test whether the visual-spatial tuning of V1 in a system that has developed using only rod input (i.e., achromats) remains stable or adapts, using population receptive field modelling. Our analyses focus on two types of plasticity. First, we ask what happens to the V1 cortex normally receiving only cone-mediated input from the foveal ‘rod-free’ zone, which shows structural abnormalities. Second, we test for changes in rod-selective visuospatial tuning across V1 (beyond the 0°-2° ‘rod-free’ fovea), in cortex which normally receives combinations of rod and cone input but develops with only rod inputs in achromats.

V1 receives retinotopically organised input, with adjacent neurons in the cortex receiving signals from adjacent neurons in the retina. Population receptive field modelling links retinotopic voxels in visual cortex to the locations in the visual field they represent (Dumoulin & Wandell, 2008). This is achieved by fitting population receptive fields estimates (pRFs; i.e., Gaussian tuning functions) to BOLD responses evoked by stimuli systematically traversing the visual field. To compare rod-mediated cortex function in achromatopsia to normally developed rod-responses under matched viewing circumstances, we used the silent substitution method (Estévez & Spekreijse, 1982) and low luminance presentation which allowed us to generate stimuli with a rod-isolating contrast (“rod-selective condition”). We also included another stimulus condition that, in a typical retina, stimulates cones and the rods equally (“non-selective condition”). This allowed us to assess how rod-mediated measures in both groups compare to more typical pRF measures driven by cones in normal sighted, and whether differences across achromats and controls remain consistent across different visual stimulus properties.

Coverage and pRF properties of the foveal projection zone

We first compared retinotopic tuning in cortex containing the rod-free foveal projection zone (atlas defined 0°-2°) across controls and achromats. In normally developed retinotopic rod maps, the cortical projection of the foveal rod-free zone, appears as a central “hole” (Figure 2B; Barton & Brewer, 2015). To test if this normally silent area fills in with reorganised retinotopic visual activity in achromats as has previously been reported (Baseler et al., 2002), we compared the percentage of cortex that contains reliable pRFs (coverage; meaning where pRFs could reliably be estimated from visually-driven signal. This is shown in Figure 2 for pRF model fits of R² > 0.1, but effects were similar at R² > 0.2 and 0.3.

Figure 2

In the rod-selective condition, the coverage in the 0°-2° region of interest was near-identical across achromats and controls (on average 24% in controls, and 25% in achromats, t₍₄₇₎ = 0.19, p=0.85, BF₁₀=0.31; see Figure 2A&B). For non-selective stimuli targeting both cones and rods, controls exhibit greater visual coverage (75% on average) compared to achromats, who are missing functional cones (43% on average) (Mann U₍₅₂₎=63, p<0.0001; unequal variance reverted to non-parametric). Achromats show enhanced coverage under non-selective compared to rod-selective conditions, likely due to increased luminance enhancing the rod system’s signal-to-noise ratio. We note that two achromats show potential filling in with a foveal activation coverage beyond the 95% predictive interval of the normal sighted rod gap. Their retinotopic profiles (see Supplement 1) suggest this may reflect individual differences in retinal foveal photoreceptor distribution profiles rather than cortical remapping. This highlights the importance of characterising variability in atypically developing populations through larger samples rather than case studies.

Crucially, data quality was matched across groups, with no differences in pRF model fits (see Supplement 2). It is therefore unlikely that any coverage inside the rod “hole” in achromats was masked by factors such as eye-movement, visual responsiveness, or head-movement. In line with this, we found no simple relationship between fixation stability and coverage (Supplement 2). In addition, it is unlikely that these findings can be explained by pRF fitting-related artefacts (Binda et al. 2013), as modelling the rod-scotoma in the pRF fitting procedure gave similar results (see methods and Supplement 3).

Together these results provide strong evidence for stability in the face of life-long deprivation, i.e., the conservation of the sensory topography. This is in line with the structural changes in this cortical area presented above, and recent conclusions by an independent group (Molz et al., 2023).

Altered rod representation in V1 - pRF size and eccentricity

We next tested for differences in spatial tuning across V1 regions that typically receive both rod and cone inputs but develop with only rod inputs in achromatopsia. Figure 2C shows that pRF sizes in achromats are systematically larger compared to controls across all eccentricities. This holds true for both the non-selective (V1 segments: 2°-3°: t₍₅₂₎=12.57, p_unc<0.001; 3°-4°: t₍₅₂₎=11.12, p_unc<0.001; 4°-5°: t₍₅₂₎=9.54, p_unc<0.001; 5°-6°: t₍₅₂₎=7.52, p_unc<0.001; 6°-7°: t₍₅₂₎=5.76, p_unc<0.001; 7°-8°: t₍₅₂₎=2.36, p_unc=0.022; all segments remain significant after FDR correction), and rod-selective stimulus conditions (V1 segments: 2°-3°: t₍₄₈₎=3.83, p_unc<0.001; 3°-4°: t₍₄₈₎=3.67, p_unc <0.001; 4°-5°: t₍₄₈₎=4.68, p_unc <0.001; 5°-6°: t₍₄₈₎=4.74, p_unc <0.001; 6°-7°: t₍₄₈₎=3.98, p_unc <0.001; 7°-8°: t₍₄₈₎=2.00, p_unc =0.051; all segments but the 7°-8° remain significant after FDR correction). In the non-selective condition, the smaller pRFs in controls are in line with the contribution the cone system’s higher spatial resolution which is not active in the achromat group. However, for the rod-selective condition where both groups activate the same class of photoreceptors, the enlarged pRFs in achromats indicate altered rod-based visual processing.

Figure 2D shows that pRF estimates in achromats have higher eccentricities across V1 compared to controls, spanning all tested V1 segments (group comparisons of the rod-selective conditions: 2°-3°: t₍₄₈₎=3.32, p_unc=0.002; 3°-4°: t₍₄₈₎=4.42, p_unc<0.001; 4°-5°: t₍₄₈₎=3.59, p_unc=0.001; 5°-6°: t₍₄₈₎=3.17, p_unc=0.003; 6°-7°: t₍₄₈₎=3.15, p_unc=0.003; 7°-8°: t₍₄₈₎=2.23, p_unc=0.030; all segments but the 7°-8° remain significant after FDR correction). However, in the non-selective condition in achromats, even though still only rods are being stimulated, estimated eccentricities are lower than those observed under the rod-selective condition. In fact, the eccentricity estimates of the non-selective condition in achromats appear similar to those of the rod-selective condition in controls. This shows that it is possible, in the achromat system, to obtain eccentricity estimates in the normal rod range. We therefore do not consider the eccentricity shifts as evidence for plasticity or stability as it is unclear whether they are a property of the tissue or reflect other factors.

Together these analyses reveal subtle differences in the rod-mediated spatial tuning of V1 cortex beyond the foveal rod-free projection zone. Notably, achromats combine information across larger areas in visual space (larger pRFs). It is not clear whether larger pRFs can be purely attributed to coarser spatial sampling at the cortical level (i.e., a degree of cortical plasticity), because this could reflect changes photoreceptor properties at the retinal level, or blurring of visual information sampling at the ocular level due to nystagmus (eye-movement). The pRF size and eccentricity differences remain unchanged when accounting for the scotoma in the pRF fitting procedure (see Supplement 3).

Connective field modelling – hierarchical sampling

We next used connective field modelling to test whether V3’s sampling along the V1 cortical surface in achromatopsia remains stable despite the drastically altered visual input from this region since birth. We test for two hypothetical types of hierarchical plasticity, shifts in the resource distribution across the V1 eccentricity map (CF centre), and shifts in the aperture of the read-out area (CF size).

Higher-order visual regions receive retinotopic information from V1, with individual neurons integrating information across multiple neighbouring V1 neurons. Connective field (CF) modelling is a functional connectivity-based extension of pRF modelling which characterises this process based on time-series correlations between brain regions. While pRF modelling is ‘stimulus referred’, estimating a voxel’s spatial tuning properties based on its response to stimulus variations, connective field modelling is ‘neural-referred’, estimating a voxel’s ‘connective field’ based on its correlations with aggregated voxels along the primary visual cortex surface. Responses throughout cortex are modelled as deriving from Gaussian shaped sampling profiles over the V1 cortical surface (connective fields), defined by a position (centre of the Gaussian) and size (Gaussian spread in mm (σ)). CF modelling can thus reveal the topographic profile of a brain region’s neural connectivity with V1, including the spatial extent of sampling. After estimating connective fields, the CF positions can be combined with the retinotopic map of V1 to ‘project’ the retinotopy of V1 onto the rest of the brain, thereby performing connectivity-derived retinotopic mapping. To create V3 connective-field-based eccentricity maps (Figure 3B), each vertex was assigned the Benson atlas-based eccentricity of its connective field centre in V1. Because CF modelling relies on inherent signal fluctuations and involves no assumptions about visual input, it is more robust to potential confounds such as poor fixation and low visual responsiveness (Tangtartharakul et al., 2023). Here we focus on connective fields between V3 and V1. We omit V2, which shares a border with V1, to avoid inaccuracies in CF estimation due to mis-delineation.

Figure 3

Shifts in V3 resource distribution along the V1 eccentricity map

When both rods and cones are stimulated in the normal sighted system, the largest proportion of V3 voxels is dedicated to sampling information from the V1 foveal projection (connective field centres in the 0°-1° degree range), with a gradual decrease in resources allocated to larger eccentricities (Figure 3A). When only rods are stimulated, foveal sampling is attenuated but the overall pattern of a gradual decrease in sampling from small to larger eccentricities remains. In achromats, the way in which V3 samples V1 is drastically different, with most V3 resources dedicated to sampling eccentricities at 1°-4° for both stimulation conditions (Figure 3B). When comparing the rod-selective condition in achromats to that of controls, we observe a shift towards greater V3 sampling from the V1 cortex encoding the (1°-4°) paracentral visual field (0°-1°: t₍₄₇₎=1.18, p_unc=0.242; 1°-2°: t₍₄₇₎=4.03, p_unc<0.001; 2°-3°: t₍₄₇₎=5.32, p_unc<0.001; 3°-4°: t₍₄₇₎=3.70, p_unc=0.001; 4°-5°: t₍₄₇₎=2.32, p_unc=0.025; 5°-6°: t₍₄₇₎=0.84, p_unc=0.405; 6°-7°: t₍₄₇₎=-0.97, p_unc=0.337; 7°-8°: t₍₄₇₎=1.93, p_unc=0.060; Segments including 1°-4° remain significant after FDR correction). It is unlikely that these differences result from differences in measurement quality, because the goodness of fit (CF R²) was comparable across the two groups (Supplement 4). This shift in sampling could reflect increased neural resources dedicated to the most fovea-adjacent information that is available to the rod-only system. Surprisingly, even though we found no evidence of visually-driven signal in the region on V1’s surface that normally takes foveal input (0°-1°) in achromats (see Figure 2), there is still a substantial portion of V3 vertices (15%) that sample from this region. This remaining foveal connectivity between V1 and V3 shows that there is also a degree of stability across retinotopic hierarchical sampling.

Shifts in the V3 read-out aperture from V1

In normal sighted, we show that connective field size in V3 increases with eccentricity, whether stimulating both rods and cones at the same time, or only rods. This replicates findings by Haak et al., (2013) and Knapen (2021). In achromatopsia, however, V3 connective field sizes remain consistently small across all eccentricities, and are significantly smaller than in the normal-sighted rod-system between 4°-6° eccentricity (comparing achromats vs control rod-selective condition: 0°-1°: t₍₄₄₎=1.07, p_unc=0.289; 1°-2°: t₍₄₅₎=-0.41, p_unc=0.684; 2°-3°: t₍₄₅₎=-0.58, p_unc=0.568; 3°-4°: t₍₄₆₎=-1.70, p_unc=0.097; 4°-5°: t₍₄₇₎=-3.55, p_unc=0.001; 5°-6°: t₍₄₇₎=-2.63, p_unc=0.011; 6°-7°: t₍₄₇₎=-0.81, p_unc=0.420; 7°-8°: t₍₄₅₎=-0.73, p_unc=0.469; degrees of freedom vary as not all participants have vertices sampling all V1 segments). These smaller connective field sizes in V3, mean that V3 voxels pool across a smaller number of pRFs along the V1 cortical surface (4C, middle).

In V1, pRF sizes in achromats were larger than in controls both for rod-selective and non-selective stimuli conditions (4C, left). Finer spatial sampling of the V1 cortex by V3 could therefore effectively compensate for this coarser input, which should result in normalisation of pRF size in V3 (4C, right). To further investigate this possibility, we plotted the ratio between pRF sizes in achromats and pRF sizes in controls, with values above 1 indicating enlarged pRFs in achromats (Figure 4B). As predicted by the smaller connective fields, the ratio of pRF sizes is closer to 1 in V3 than in V1, indicating a normalisation of pRF size in the higher-order region relative to the input stage (Figure 4D). This is consistently the case in both stimulus conditions. This finer hierarchical sampling could reflect a higher-order compensation for coarser visual inputs in achromats.

Figure 4

Discussion

To test the effects of long-term altered visual input on neural plasticity across the cortical hierarchy, we compared a unique group of individuals who have developed only with rod photoreceptor inputs (achromatopsia) to normally sighted individuals who spend most of their waking time relying on cone-dominated vision. As normal-sighted can also experience rod-only vision under specific circumstances, this allowed us to dissociate the effects of life-long altered experience from transient altered sensory input. By integrating markers of plasticity and stability across structural, representational, and connectivity analyses we reveal system-level dynamics of stability and adaptations to altered visual experience. The results offer a new perspective on neural plasticity and stability, recasting these processes as two interacting forces acting simultaneously on the cortex rather, than two separate neural responses to change.

Specifically, in achromatopsia we found thickening of the cortex around the foveal projection zone, which receives only cone-photoreceptor inputs and therefore would be signal-deprived from birth in this group. In contrast to early reports (Baseler et al., 2002), we found no evidence that visually-induced signalling inside this region had increased compared to the typically developed rod system. We also observed subtle functional changes in V1 retinotopic tuning including increased population receptive field sizes, which may reflect retinal or post-retinal changes. These data reveal that visual spatial representations in V1 are highly stable, integrating and resolving ambiguity of earlier findings (Baseler et al., 2002; Mckyton et al., 2021; Molz et al., 2022, 2023). This high stability at the cortical input-stage contrasted with reorganisation of the read-out of these inputs by higher-order visual cortex. Notably, more cortical resources in V3 were dedicated to sampling input from the near-fovea, and V3 neuron populations aggregated information across smaller areas along the V1 cortical surface. There was also evidence for stability of sampling in achromatopsia, because V3 retained connectivity with the foveal confluence despite the life-long scotoma.

Together these findings show that while primary sensory areas such as V1 exhibit considerable stability, higher-order areas such as V3 display marked plasticity that may optimise use of the available information. Both increased sampling of the most fovea-adjacent information available to rods, and finer-grained sampling of coarser inputs could effectively counteract the limitations of developing with a fully rod-based system. These findings add a new dimension to the ongoing debate on neural plasticity versus stability by proposing a hierarchical network perspective where both processes interact across processing stages to maintain the integrity of foundational architecture, whilst optimizing the processing of available input.

Control measures and analyses presented above rule out several alternative explanations for these results. Model fits were similar for both patients and controls, so group differences are unlikely to reflect data quality discrepancies. Further, eye-tracking data, available for a subset of participants, revealed no correlation between stable fixation and enhanced foveal activity in achromatopsia (Supplement 2), suggesting that ’filling in’ of the deprived V1 cortex was not simply obscured by nystagmus. This is further supported by converging structural evidence of cortical thickening of this region, a typical outcome of congenital visual deprivation (Jiang et al., 2009; Park et al., 2009). Changes in V3 sampling of V1 are also unlikely to be a by-product of eye-movement, as connective field modelling is locked to V1 inputs rather than screen locations, so insensitive to eye movement (Tangtartharakul et al., 2023).

There are convincing arguments for the importance of a system that maintains highly stable sensor-to-cortex connections at the first cortical input stage. Re-routing of cortical projections in response to changing visual signals could lead to “coding catastrophe”, because left unchecked it could leave visual processing areas throughout the brain with mislabelled retinotopic signals, corresponding to different regions of visual space than they would normally represent (Dumoulin & Knapen, 2018; Haak et al., 2015; Wandell & Smirnakis, 2010). Higher-order regions aggregate more inputs across visual space and other sensory domains, and display greater invariance to low-level stimulus properties (Lehky & Tanaka, 2016). Therefore, there is theoretically broader scope for experience-dependent reweighting of inputs (Makin & Krakauer, 2023), to optimise use of inputs that are still available, more reliable, or more relevant in the impaired system.

We propose a new hierarchical model of neural adaptation that hypothetically could act to improve the limitations of a rod-only visual system (i.e., lacking foveal input and reduced acuity). However, to confirm that the adjustments we observe are compensatory rather than mere by-products of circuitry responses to visual deprivation, their functional benefits will need to be demonstrated. The adaptations might suggest that individuals with achromatopsia exhibit superior rod-based vision compared to typical systems. One study supporting this, reported improved rod-mediated biological motion processing in some achromatopsia patients, suggesting an enhanced rod system capacity (Burton et al., 2016). Nonetheless, most studies typically find achromatopsia patients’ rod vision to be similar or inferior, potentially because retinal structure can also be impacted by the disease (Burton et al., 2016; Hess & Nordby, 1986). It is likely that the functional outcome is shaped by a mixture of impairments and compensatory mechanisms at various levels, presenting an ideal case for examining such dynamics.

In light of breakthrough treatments set up to restore cone function in human patients with achromatopsia (Farahbakhsh et al., 2022; Fischer et al., 2020; Michaelides et al., 2023; Reichel et al., 2022), the new plasticity dynamics revealed by our work have significant theoretical and practical implications. The possible dormancy of retinotopic maps in cortical input stages as well as some preserved connectivity between V1 foveal areas and V3 bodes well for these gene therapies, each indicating remaining infrastructure to process rescued cone signals. However, our findings also show that the interpretation of signals from V1 by higher-order networks is modified, and it remains to be seen how this effects the translation of new sensory inputs into functional vision. Adaptations in neuronal tuning properties vary in how stably ingrained they become in neural architecture. They can temporarily shift in response to task-demands (transient modulation) (e.g., Klein et al., 2018), or unfold over longer periods with extensive experience that can be acquired or reversed at any age (perceptual learning) (Lu & Dosher, 2022). They can also involve more profound reweighting of neural properties that do not easily revert after critical developmental windows (developmental divergence) (Mitchell & Maurer, 2022). A previous study found that connective fields stay stable with central field loss in adulthood (Haak et al., 2016), raising the possibility of early windows of plasticity for these adaptations. Understanding whether the adaptations in achromatopsia are within the normal range of flexibility for a typically developed system or represent a fundamental reorganisation will be of high relevance to treatment development, as it impacts the nature and extent of neural plasticity that can be harnessed for therapeutic recovery at different ages.

With over 50 interventional studies for inherited retinal diseases initiated, the technology to treat currently incurable heritable retinal diseases by restoring ocular cellular structure and function is here. The next big hurdle will be to understand neural plasticity constraints on integrating new signals from the eye into the existing neural infrastructure. Our research provides a novel framework for understanding the interactions of neural stability and plasticity as the human brain adapts to altered information, which is now of crucial importance for harnessing the full potential of emerging visual restoration therapies.

Methods

Participants

All participants met MRI safety inclusion criteria and had no other known neurological disorders. Patients with achromatopsia had genetically confirmed disease (CNGA3 or CNGB3 bi-allelic disease-causing variants). The main dataset was collected at 1.5T. Additional structural data was collected at 3T.

Main 1.5T dataset

We obtained functional and structural 1.5T data from fourteen individuals with achromatopsia (age range = 8-27; mean±std = 13.7±5.7), and forty normally sighted controls (age range=6-34; mean±std=19.8±8.3). Additional structural 1.5T data was collected from 8 achromats (age range=8-35; mean±std=19.75±10.6) and 50 normally sighted controls (age range=6-29; mean±std=12.5±5.9).

Additional 3T dataset (previously published in Anderson et al. 2024)

We obtained additional structural scans from ten normally sighted controls (age range=20-30; mean±std=23.8±3.7) and 9 individuals with achromatopsia (age range=16-33; mean±std=24.2±6.4).

For 4 normally-sighted individuals the rod-only pRF map was not collected. One achromat was excluded from functional analysis due to large head-movements. Data collection had ethics approval from the national ethics committee for individuals with achromatopsia (REC reference: 12/LO/1196, IRAS code: 106506, REC reference 11/LO/1229) and the UCL Research Ethics Committee for normal sighted control participants (4846/001).

MRI acquisition

Structural MRI

All 1.5T participants were scanned in a Siemens Avanto 1.5T MRI scanner with a 30-channel coil (a 32-channel coil customised to remove view obstructions). A high-resolution structural scan was acquired using a T1-weighted 3D MPRAGE (1 mm³ voxel size, Bandwidth = 190 Hz/Px, 176 partitions, partition TR = 2730, TR = 8.4ms, TE = 3.57, effective T1 = 1000 ms, flip angle = 7 degrees).

All 3T participants were scanned in a Siemens 3T TIM-Trio scanner using a 32-channel head coil. A high-resolution structural scan was acquired using a T1-weighted 3D MDEFT, 1mm isotropic voxels, 176 sagittal slices, 256x240 matrix, TE 2.48ms, TR 7.92ms, TI 910ms.

Functional MRI

Functional T2*-weighted multiband 2D echo-planar images were collected using a multi-band sequence with 2.3 mm isotropic voxels (TR = 1000 ms, TE = 55 ms, volumes = 348, flip angle = 75°, MB acceleration factor = 4, Bandwidth = 1628 Hz/Px).

pRF stimuli and apparatus

Stimuli were presented on an MR-compatible LCD display (BOLDscreen 24, Cambridge Research Systems Ltd., UK; 51 x 32 cm; 1920 x 1200 pixels) viewed through a mirror in the scanner at 105 cm distance. Participants were lying supine in the scanner, with fixation stability recorded where possible via a mirror, with an Eyelink 1000 Plus (SR Research, Ottawa, ON) at the back of bore. Built-in gaze calibration was not possible in achromats due to nystagmus so the camera was calibrated in advance on a healthy eye. Before every other run, patients fixated on a 5-point custom calibration which allowed us to calibrate gaze measures post-hoc (Tailor et al., 2021) Stimuli were presented using MATLAB (R2016b, MathWorks, Natick, MA, USA), via the Psychophysics Toolbox 3 (Kleiner et al., 2007).

pRF mapping stimuli comprised a simultaneous rotating ring and contracting/expanding wedge embedded within a contrast-reversing checkerboard with 2Hz reversal rate (Figure 1B). Per run, the ring expanded/contracted for 6 cycles (48 secs/cycle) with logarithmic eccentricity scaling (van Dijk et al., 2016). The wedge (20° angle) rotated clockwise/anticlockwise for 8 cycles (36 sec/cycle). 20-second fixation baselines were embedded at the start, at the mid-point, and end of the run (total run duration 348 secs). The stimuli covered a maximum eccentricity of 8.6°, and moved to a new position each 1-second TR. The stimulus was overlaid with a small white central fixation dot (0.2° VA radius) and a black radial grid to encourage stable fixation. To keep participants engaged, they played a rewarded “kitten rescue mission game” in which they pressed a button when detecting a fixation target luminance change.

The rod-selective (stimulating only rods) and non-selective stimuli (stimulating rods and cones equally) differed in both the colour pair (chromatic pair) used to make up the checkerboard and luminance levels. For an in depth description of stimulus generation and validation see Farahbakhsh et al., 2022. In brief, a silent substitution approach was used to obtain a rod-selective colour pair (Estévez & Spekreijse, 1982). We generated a single chromatic pair that kept L and M cone contrast constant (i.e., invisible to L and M cones), whilst inducing a contrast response in rods. Since it is only possible to simultaneously silence two photoreceptor types with 3 colour channels, the resulting stimulus was not controlled for S-cone contribution. To eradicate the S-cone response, we presented these stimuli at very low light level (max luminance: 0.02 cd/m²). The non-selective stimuli chromatic pair increased rod and L and M cone contrasts equivalently. This stimulus was presented at mesopic luminance levels (max 0.5 cd/m2).

Prior to the rods only condition, which always appeared before the non-selective condition participants laid in the scanner for a 15-minute dark adaptation. A group of 10 normal sighted control adults had longer period of dark adaptation (45 minutes) test any impact on the resultant pRF maps. However, since we did not find any differences between the two groups they were combined and treated as a single control group. We note that excluding these 10 participants did not change any of the reported results.

MRI pre-processing

Data were preprocessed using fMRIPrep 21.0.2 (Esteban, Blair, et al., 2018; Esteban, Markiewicz, et al., 2018; RRID:SCR_016216), which is based on Nipype 1.6.1 (K. Gorgolewski et al., 2011; K. J. Gorgolewski et al., 2018; RRID:SCR_002502).

fmriprep structural data pre-processing

The T1-weighted (T1w) image was corrected for intensity non-uniformity (INU) with N4BiasFieldCorrection (Tustison et al., 2010), distributed with ANTs 2.3.3 (Avants et al., 2008; RRID:SCR_004757), and used as T1w-reference throughout the workflow. The T1w-reference was then skull-stripped with a Nipype implementation of the antsBrainExtraction.sh workflow (from ANTs), using OASIS30ANTs as target template. Brain tissue segmentation of cerebrospinal fluid (CSF), white-matter (WM) and gray-matter (GM) was performed on the brain-extracted T1w using fast (FSL 6.0.5.1:57b01774, RRID:SCR_002823, Zhang et al., 2001). Brain surfaces were reconstructed using recon-all (FreeSurfer 6.0.1, RRID:SCR_001847, Dale et al., 1999), and the brain mask estimated previously was refined with a custom variation of the method to reconcile ANTs-derived and FreeSurfer-derived segmentations of the cortical gray-matter of Mindboggle (RRID:SCR_002438, A. Klein et al., 2017). Volume-based spatial normalization to one standard space (MNI152NLin2009cAsym) was performed through nonlinear registration with antsRegistration (ANTs 2.3.3), using brain-extracted versions of both T1w reference and the T1w template. The following template was selected for spatial normalization: ICBM 152 Nonlinear Asymmetrical template version 2009c [Fonov et al., 2009, RRID:SCR_008796; TemplateFlow ID: MNI152NLin2009cAsym]. Many internal operations of fMRIPrep use Nilearn 0.8.1 (Abraham et al., 2014, RRID:SCR_001362), mostly within the functional processing workflow.

fmriprep functional data pre-processing

For each BOLD run found per subject (across all tasks), the following preprocessing was performed. First, a reference volume and its skull-stripped version were generated by aligning and averaging 1 single-band references (SBRefs). Head-motion parameters with respect to the BOLD reference (transformation matrices, and six corresponding rotation and translation parameters) are estimated before any spatiotemporal filtering using mcflirt (FSL 6.0.5.1:57b01774, Jenkinson et al., 2002). The BOLD time-series were resampled onto their original, native space by applying the transforms to correct for head-motion. These resampled BOLD time-series will be referred to as preprocessed BOLD in original space, or just preprocessed BOLD. The BOLD reference was then co-registered to the T1w reference using bbregister (FreeSurfer) which implements boundary-based registration (Greve & Fischl, 2009). Co-registration was configured with six degrees of freedom. First, a reference volume and its skull-stripped version were generated using a custom methodology of fMRIPrep. Several confounding time-series were calculated based on the preprocessed BOLD: framewise displacement (FD), DVARS and three region-wise global signals. FD was computed using two formulations following Power (absolute sum of relative motions, Power et al., 2014) and Jenkinson (relative root mean square displacement between affines, Jenkinson et al., 2002). FD and DVARS are calculated for each functional run, both using their implementations in Nipype (following the definitions by Power et al., 2014). The three global signals are extracted within the CSF, the WM, and the whole-brain masks.

Additionally, a set of physiological regressors were extracted to allow for component-based noise correction (CompCor, Behzadi et al., 2007). Principal components are estimated after high-pass filtering the preprocessed BOLD time-series (using a discrete cosine filter with 128s cut-off) for the two CompCor variants: temporal (tCompCor) and anatomical (aCompCor). tCompCor components are then calculated from the top 2% variable voxels within the brain mask. For aCompCor, three probabilistic masks (CSF, WM and combined CSF+WM) are generated in anatomical space. The implementation differs from that of Behzadi et al. in that instead of eroding the masks by 2 pixels on BOLD space, the aCompCor masks are subtracted a mask of pixels that likely contain a volume fraction of GM. This mask is obtained by dilating a GM mask extracted from the FreeSurfer’s aseg segmentation, and it ensures components are not extracted from voxels containing a minimal fraction of GM. Finally, these masks are resampled into BOLD space and binarized by thresholding at 0.99 (as in the original implementation). Components are also calculated separately within the WM and CSF masks. For each CompCor decomposition, the k components with the largest singular values are retained, such that the retained components’ time series are sufficient to explain 50 percent of variance across the nuisance mask (CSF, WM, combined, or temporal). The remaining components are dropped from consideration. The head-motion estimates calculated in the correction step were also placed within the corresponding confounds file. The confound time series derived from head motion estimates and global signals were expanded with the inclusion of temporal derivatives and quadratic terms for each (Satterthwaite et al. 2013). Frames that exceeded a threshold of 0.9 mm FD or 1.5 standardised DVARS were annotated as motion outliers. The BOLD time-series were resampled onto the following surfaces (FreeSurfer reconstruction nomenclature): fsnative. All resamplings can be performed with a single interpolation step by composing all the pertinent transformations (i.e., head-motion transform matrices, susceptibility distortion correction when available, and co-registrations to anatomical and output spaces). Surface resampling was performed using mri_vol2surf (FreeSurfer).

Additional pre-processing

The first 4 volumes of each run were discarded for signal stabilization. Bold data in fsnative space was smoothed along the cortical surface (Gaussian kernel FWHM = 3 mm). For each run, a denoising procedure was further applied on the data by inputting nuisance regressors into a PCA and regressing out the top 5 components. The nuisance regressors include: Six head motion regressors, FD and DVARS, CSF mean signal, white matter mean signal, the first five anatomical CompCor regressors (csf+wm combined), and the first five temporal CompCor regressors. Then data and regressors were temporally filtered with a Savitzky-Golay filter (3rd order, 347 seconds in length) and zscored. The time series from the two identical task runs were then averaged to increase signal to noise ratio.

MRI analysis

Region of interests

An anatomically defined retinotopy atlas was applied to the inflated cortical surface for each participant in order to delineate the borders of V1, V2 and V3, and also each degree of eccentricity up to 8° in V1 (Benson & Winawer, 2018). See Figure 1C for an example of an eccentricity plot for one participant.

Structural data analysis

The ‘mris_anatomical_stats’ command was used to extract surface area and cortical thickness properties for each participant, per hemisphere, for each region of interest. Overall cortical thickness per participant was computed by weighting the cortical thickness of a hemisphere by the number of vertices per hemisphere.

pRF analysis

pRF modelling was carried out with a custom Python package, prfpy (Aqil & Knapen, 2023). To model population receptive fields, we used a isometric bivariate Gaussian model (Dumoulin & Wandell, 2008), with x,y representing the preferred retinotopic position, spatial standard deviation sigma (σ) representing pRF size as well as the baseline of the timeseries and the signal amplitude beta (β). At every cortical location, the parameters of each model were optimized by minimizing the residual sum of squares between the model prediction and the measured BOLD signal. Fitting began with a grid search for the best model parameters (pRF position and size) followed by an iterative optimization. In each participant, in each condition, vertices with an R² < 0.1 (variance explained by the pRF model), a negative beta, or a centre outside the stimulated radius (8.6°) were excluded from further analysis.

In addition, we ran a second model fit using a stimulus aperture file that incorporated a black oval mask (1.25°x0.9°) over the area of the central rod-free scotoma (Curcio et al., 1990). Previous research suggests that biases in pRF measures such as shifts in eccentricities and size can occur when a scotoma is present and not explicitly modelled in the stimulus apertures (Binda et al., 2013). Results from this analysis are reported in Supplement 3.

Connective field modelling

Connective field modelling was carried out with a dedicated Python package, prfpy (Aqil & Knapen, 2023) based on the method described in Knapen (2021). In brief, gaussian connective field profiles on the surface are defined for each vertex on the cortical manifold, with a centre location of the connective field (v-centre), and σ (sigma) the Gaussian spread in millimeters on the cortical surface.

To generate spatial Gaussian fields along the V1 surface, the diffusion of heat along the cortical mesh was simulated (Crane et al., 2013) and used to infer geodesic distances between all V1 vertices (V1 definition was based on the Benson atlas, see Regions of Interest). As the two hemispheres are two separate surface meshes, this distance matrix was calculated for each hemisphere separately. Grid-fit candidate CF sizes were: [0.5, 1, 2, 3, 4, 5, 7, 10, 15, 20] mm σ for the Gaussian CF. The predicted time course for each of the resulting CF models was generated by taking the dot product between the CF’s vertex profile and the vertex by time-courses matrix in V1. The resulting CF model time courses were z scored and correlated with the time courses throughout the brain. For each location in the brain, the connective field resulting in the highest squared correlation was selected. Only values within an R² > 0.1, beta > 0, and a centre within 8° of eccentricity (matching the radius of the stimuli) were included in further analyses. To create connective-field-based eccentricity maps, each vertex was assigned the Benson atlas-based eccentricities (Benson & Winawer, 2018) of its CF V1 centre.

For visualisation of group maps all individual maps were transformed to fsaverage space. Plotted values indicate the medians of each group and condition per vertex. A vertex was left empty if no significant data was available for more than 30% of individuals within the group.

Statistical analysis

Jamovi R package (jmv) was used for Anovas and t-tests (The jamovi project, 2022). Jasp (jsq) R package for Bayesian t-tests (JASP team, 2020). Where sphericity assumptions were not met a greenhouse-geisser correction was applied. Uncorrected p-values are reported as well as FDR multiple comparison corrections.

Acknowledgements

We thank Maria Molina Sanchez, Hunter Schone and Matan Mazor for providing comments on an earlier version of this manuscript.

Funding

The research was supported by grants from the National Institute for Health Research (NIHR) Biomedical Research Centre (BRC) at Moorfields Eye Hospital NHS Foundation Trust and UCL Institute of Ophthalmology, the Economic and Social Research Council (ESRC) of the UKRI (#ES/N000838/1); MeiraGTx, Retina UK, Moorfields Eye Charity (R160035A, R190029A, R180004A), Foundation Fighting Blindness (USA) and The Wellcome Trust (099173/Z/12/Z, 100227), Ardalan Family Scholarship, the Persia Educational Foundation Maryam Mirzakhani Scholarship; and the Sir Richard Stapley Educational Trust (#313812). The views expressed are those of the authors and not necessarily those of the NHS, the NIHR, or the Department of Health.

Competing interests

Scanner time was funded by MeiraGTx via a research sponsorship awarded to UCL. M.M. declares financial interest in MeiraGTx.