Cortical quantity representations of visual numerosity and timing overlap increasingly into superior cortices but remain distinct

Many sensory brain areas are organized as topographic maps where neural response preferences change gradually across the cortical surface. Within association cortices, 7-Tesla fMRI and neural model-based analyses have also revealed many topographic maps for quantities like numerosity and event timing, often in similar locations. Numerical and temporal quantity estimations also show behavioral similarities and even interactions. For example, the duration of high-numerosity displays is perceived as longer than that of low-numerosity displays. Such interactions are often ascribed to a generalized magnitude system with shared neural responses across quantities. Anterior quantity responses are more closely linked to behavior. Here, we investigate whether common quantity representations hierarchically emerge by asking whether numerosity and timing maps become increasingly closely related in their overlap, response preferences

Here we investigate the relationships between neural responses to different quantities that may underlie these behavioral interactions.At the neural level, similar networks of brain areas, including areas around the intraparietal sulcus (IPS), change their response pattern with different values of each quantity (Borghesani et al., 2019;Eger et al., 2009;Harvey et al., 2020;Harvey and Dumoulin, 2017a;Hayashi et al., 2018).Indeed, intraparietal activity seems necessary for discrimination within each quantity (Alexander et al., 2005;Cohen Kadosh et al., 2007;Dormal et al., 2012;Lecce et al., 2015;Riemer et al., 2016).Recently, ultra-high-field (7T) fMRI and neural model-based analyses have revealed similarly located networks of cortical maps with a topographic progression of preferred quantities (Harvey et al., 2015), including the numerosity of visual objects (Harvey et al., 2013;Harvey and Dumoulin, 2017a) and the timing (duration and period) of visual events (Harvey et al., 2020) (also see Protopapa et al., 2019).Evidence that such quantity-selective responses are linked to perception of their respective quantity is accumulating (reviewed by Tsouli et al., 2022).
Similarities and interactions between perception of different quantities, together with similar brain areas responding to tasks involving different quantities, have led to the prominent 'A Theory Of Magnitude'.This suggests that different quantities are transformed onto a single common metric of quantity that is used to guide motor responses to all quantities (Walsh, 2003).Interpreted most literally, this predicts that responses to different quantities would be just as closely related as repeated measures of responses to the same quantity.However, an alternative explanation is that different quantities produce (at least partially) distinct responses within the same brain area, allowing interactions between the neural populations responding to each quantity (Hofstetter et al., 2021;Tsouli et al., 2022).
Potentially, responses to different quantities might become increasingly related through a hierarchical transformation.Indeed, numerosity and timing are initially derived from different sensory information at different processing stages (Hendrikx et al., 2022;Paul et al., 2022), but have responses in similar areas later in their respective hierarchies.Furthermore, within each quantity we also see hierarchical transformations, such that from posterior to anterior brain areas, neural response preferences for both quantities become more specific (Harvey et al., 2020), more tightly linked to behavior (Ramirez-Cardenas et al., 2016;Viswanathan and Nieder, 2015), and more separated from non-quantity sensory properties like visual position (Hendrikx et al., 2022;Paul et al., 2022).These properties raise possibility that numerosity and timing may initially be estimated through separate processing of sensory inputs and then brought together in premotor areas to support the planning of common motor responses, as A Theory Of Magnitude predicts.
The relationships between neural responses to different quantities is complicated by the fact that some quantities often co-occur in everyday situations.For example, the haptic and visual numerosity of a single group of items held in the hand is the same.This alone could underlie perceptual interactions (due to an expectation of a relationship between these quantity pairs) and the similarity of neural responses (due to neural populations responding to the same haptic and visual numerosity usually being activated together, potentially strengthening their interactions).Indeed, neural numerosity preferences are reliably correlated between visual and haptic modalities, though less strongly correlated than repeated measures of numerosity preferences within the visual or within the haptic modality (Hofstetter et al., 2021).
Other pairs of quantities, like visual numerosity and event timing, are not related so strongly (if at all) in natural stimuli.Even in such cases, perceptual interactions are still found (Arrighi et al., 2014;Cappelletti et al., 2009;Lu et al., 2009;Oliveri et al., 2008;Tsouli et al., 2019) but could depend on related actions being taken to the different quantities, particularly in laboratory situations where participants are often asked to respond to comparisons within each quantity using the same buttons (Dramkin et al., 2022;Picon et al., 2019).Likewise, neither numerosity or timing is strongly related to the spatial location of stimuli, but there are still some interactions between both quantities and stimulus location (De Tommaso and Prpic, 2020;Dehaene et al., 1993;Ishihara et al., 2008), which have also been proposed to reflect response selection biases.Similar brain networks are activated for many quantities and these areas also contain spatial representations (Bueti and Walsh, 2009).In these areas, neural spatial preferences and quantity preferences appear unrelated (Harvey et al., 2015;Harvey and Dumoulin, 2017a;Hendrikx et al., 2022;Paul et al., 2022;Viswanathan and Nieder, 2020).It remains unclear whether the same neural populations within these areas respond to multiple quantities and whether these similar response locations reflect responses to the same comparison task (Harvey et al., 2017), rather than shared responses to different stimulus quantities.Examining the relationship between neural responses to such seemingly unrelated quantities, in the absence of related quantity comparisons tasks, therefore allows testing directly whether different quantities are transformed onto a common metric.
Here we ask whether numerosity and event timing are being hierarchically transformed into a common representation of quantities that might underlie their behavioral interactions.We address this question by reanalyzing 7T fMRI data that was acquired in two previous studies mapping neural responses to numerosity and event timing (Harvey et al., 2020;Harvey and Dumoulin, 2017a).The numerosity stimuli comprised a number of circles presented on a screen (Harvey et al., 2013;Harvey and Dumoulin, 2017a).The timing stimuli comprised repetitive visual events which gradually varied in event duration and/or period (Harvey et al., 2020).To quantify the relationship between responses to these quantities, we compared the overlap between the numerosity and timing maps throughout their respective hierarchies, as well as relationships between neural numerosity and timing response preferences and their topographic progressions across the cortical surface.If different quantities were increasingly transformed onto a common metric, we would expect the overlap between cortical maps for different quantities to hierarchically increase until it matches the overlap between repeated measures of maps for the same quantity.This transformation would likewise predict that the correlation between neural quantity preferences and the similarity between direction of topographic map progressions would increase to match that of repeated measures.We found that the cortical overlap between numerosity and timing maps increased inferior-to-superior following the brain's hierarchical transformation of quantity preferences, but we did not find consistent or increasing relationships between numerosity and timing preferences.

Methods
The current study is a reanalysis of existing data from previous publications (Harvey et al., 2020;Harvey and Dumoulin, 2017a).The Methods described before the section 'Overlap of numerosity and timing maps' are shared with parts of these studies.

Participants
We collected data from eight right-handed participants (1 female, aged 25 to 35) with normal or corrected to normal vision, as described previously (Harvey et al., 2020;Hendrikx et al., 2022).This sample was a convenience sample, including two co-authors aware of the study's objectives.Other participants were financially compensated for their time and travel costs.All participants provided written informed consent.All experimental procedures were cleared by the ethics committee of University Medical Center Utrecht (protocols: 09/350 and 17/414).

Stimuli
All visual stimuli were generated using MATLAB and Psychtoolbox (Brainard, 1997;Kleiner et al., 2007;Pelli, 1997).Visual stimuli were presented on a 27.0 × 9.5 cm screen inside the MRI bore (resolution E. Hendrikx et al. 1600 × 538 pixels) at 41 cm from the participant's eyes.Participants were asked to fixate at the center of a red fixation cross that crossed the whole display on a gray background.
Quantities were presented on a gray background within 0.75 • (radius) around the center of a large, red fixation cross, as described previously (Harvey et al., 2015(Harvey et al., , 2020;;Harvey and Dumoulin, 2017a;Hendrikx et al., 2022;Paul et al., 2022).For both quantities, participants were asked to indicated with a button press when the circles were white, rather than the usual black.This happened pseudo-randomly, equally frequently for all numerosities or timings (~10 % of the trials).No quantity judgements were required.
Before scanning, participants underwent a short training in the particular quantity discrimination tasks to minimize learning or habituation effects during scanning and encourage attention to the stimulus.

Numerosity stimuli
Each numerosity pattern was displayed for 300 ms, followed by a blank gray background for 400 ms, followed by a different numerosity pattern.Six different patterns with the same numerosity were presented in a row, such that two TRs (each 2100 ms) contained one specific numerosity.Numerosities 1 through 7 were first presented in ascending order.This was followed by 16.8 s of numerosity 20, which allowed us to Fig. 1.Neural response model fitting procedure.(A) Conceptual summary.Each model began with a large set of candidate parametric response functions with different parameters.These were evaluated at every presented quantity state at the time it was presented to give a candidate set of neural response predictions: the neural response time course predicted by each candidate parameter set, given the stimulus time course.These were each convolved with a hemodynamic response function to give a candidate set of fMRI response predictions.Each prediction's correlation with the response time course of every fMRI voxel was then calculated to find the candidate response function parameters that best predicted the response measured at every voxel.(B) Numerosity response model fitting procedure, showing the response at a single fMRI voxel and its best fitting candidate response function.Here, the numerosities were each presented for the same amounts of time and the response function was one-dimensional.(C) Timing response model fitting procedure, showing the response at a single fMRI voxel and its best fitting candidate response function.Here, several stimulus configurations (see Methods) were used to sample the timings in a two-dimensional space.These configurations each produce a one-dimensional slice through the response function.In each configuration, events with different timings were each presented at different frequencies.To account for this, we evaluated the response function at the offset of every event to predict the neural response time course.Because we do not expect neural response amplitude to increase linearly with event frequency, we add a further parameter (frequency exponent) to capture this nonlinearity.When the resulting neural response time courses for each stimulus configuration are convolved with the hemodynamic response function, this produces four predictions for the response to the four stimulus configurations.For each voxel, we find the parameters of the response function and frequency exponent that minimize the different between all the predictions and responses across all stimulus configurations.E. Hendrikx et al. distinguish numerosity-selective responses with very small tuning widths (which would respond briefly and with low amplitude to numerosity 1 through 7) from those with very large tuning widths (which respond continuously and with high amplitude to numerosity 1 through 7) in fMRI responses (Dumoulin and Wandell, 2008;Harvey et al., 2013).Numerosities 7 through 1 were then presented in descending order, to counterbalance for adaptation effects of preceding numerosities, followed by another 16.8 s of numerosity 20 (Fig. 1B).This stimulus cycle was repeated four times, resulting in scanning runs of 369.6 s in total, typically for 8 runs per scan session.
For participants 1 and 2 (of 8), we systematically varied nonnumerical stimulus properties in four configurations to distinguish responses to numerosity from responses to other low-level stimulus properties (Supplementary Fig. 1; Supplementary Movie 1).In the first configuration, constant area, the summed area of all circles, and thus the mean luminance of the display, remained constant across numerosities.In the second configuration, constant object size, the size of the individual circles remained constant across numerosities.In the third configuration, constant perimeter, the summed perimeters of all circles, and thus length of all edges in the display, remained constant across.In the fourth configuration, high density, circles were presented close to each other, within a 0.375 • radius area that was placed randomly on every display to still fall entirely within 0.75 • of fixation.In other configurations the positions of the circles were spread evenly and pseudorandomly.Having found very similar responses in all these configurations, we only used the constant area and constant size configurations for the remaining six participants (3 to 8).Responses to these configurations were scanned in different sessions.

Timing mapping stimuli
Timing mapping stimuli varied in both event duration (the time between appearance and disappearance of the circle) and event period (the time between the appearance of consecutive circles, i.e., 1/frequency).Each event timing was repeated for 2100 ms, such that one fMRI volume acquisition (TR, 2100 ms) contained one specific combination of duration and period.However, the change in event timing did not always exactly coincide with a new TR since 2100 ms is not an exact multiple of all events' periods.The maximum difference between the new event timing and new TR was only 300 ms, and changes in event duration and period were in 50 ms steps, so this deviation was not perceptible.The actual presented timings are used in our modeling procedure.
The relation between event duration and period was systematically varied in four configurations to distinguish responses to event timing from responses to other low-level stimulus properties (Supplementary Movie 2, Fig. 1C).In the constant luminance configuration, event duration and period were always equal i.e., there was always a circle on the screen (Fig. 1C, black points).Duration and period first increased in 50 ms steps from 50 ms to 1000 ms, followed by a duration of 2000 ms and period of 2100 ms for 16.8 s that again allowed us to distinguish timing-selective responses with very small and very large tuning widths.Duration and period then decreased from 1000 ms to 50 ms, again counterbalancing for adaptation effects of preceding timings, again followed by a duration of 2000 ms and period of 2100 ms for 16.8 s.
In the constant duration configuration, event duration was fixed at 50 ms while event period varied (Fig. 1C, blue points).Event period first increased in 50 ms steps from 50 ms to 1000 ms, followed by a period of 2100 ms for 16.8 s.Event period then decreased from 1000 ms to 50 ms, again followed by a period of 2100 ms for 16.8 s.
In the constant period configuration, event period was fixed at 1000 ms while event duration varied (Fig. 1C, red points).Event duration first increased in 50 ms steps from 50 ms to 1000 ms, followed by a duration of 2000 ms and period of 2100 ms for 16.8 s.Event duration then decreased from 1000 ms to 50 ms, again followed by a duration of 2000 ms and period of 2100 ms for 16.8 s.
In the gaps configuration, previously unsampled combinations were presented (Fig. 1C, green points).Here each timing progression was shorter, sampling 10 timing states rather than 20, again changing in ms steps.We first presented events increasing duration from 50 ms to 500 ms while decreasing period from 950 ms to 500 ms, followed by 6.3 s with 50 ms duration and 2100 ms period, followed by events increasing duration from 50 ms to 500 ms while increasing period from 550 ms to 1000 ms, followed by 6.3 s with 50 ms duration and 2100 ms period, followed by events decreasing duration from 500 ms to 50 ms while increasing period from 500 ms to 950 ms, followed by 6.3 s with 50 ms duration and 2100 ms period, followed by events decreasing duration from 500 ms to 50 ms while decreasing period from 1000 ms to 550 ms, followed by 14.7 s with 50 ms duration and 2100 ms period.These four stimulus configurations were presented in the same scan run in all possible orders, resulting in 24 scan runs per participant, acquired in four sessions.Each scan run lasted 470.4 s in total.

fMRI data collection and pre-processing
All scans were acquired using a 7T Philips Achieva scanner, as previously described (Harvey et al., 2015(Harvey et al., , 2020;;Harvey and Dumoulin, 2017a;Hendrikx et al., 2022;Paul et al., 2022).T1-weighted anatomical scans were segmented automatically using Freesurfer.This segmentation was manually edited to reduce segmentation errors using ITK-SNAP.The T2*-weighted functional scans were acquired using a 32-channel head coil at a resolution of 1.77 × 1.77 × 1.75 mm, with interleaved slices of 128 × 128 voxels, resulting in a field of view of × 227 × 72 mm.This omitted the anterior temporal and frontal lobes where 7T fMRI has a low response amplitude and substantial spatial distortions.The repetition time (TR) was 2100 ms, echo time (TE) 25 ms, and flip angle 70 • .We used a single shot gradient echo sequence with SENSE acceleration factor 3.0 and anterior-posterior encoding.We used a 3rd-order image-based B0 shim of the functional scan's field of view (in-house IDL software, v6.3, RSI, Boulder, CO, USA).With these functional scans, we also acquired a T1-weighted anatomical image with the same resolution, position and orientation as the functional data.For numerosity mapping, we used 8 functional runs of each 176 TRs (369.6 s) for each stimulus configuration.For timing mapping, we used 24 runs of 224 TRs (470.4 s).Between the timing scans, we also acquired top-up scans with the opposite phase-encoding direction in the gradient encoding direction, to facilitate correction for image distortion during co-registration of anatomical and functional (Andersson et al., 2003) data, as previously described (Harvey et al., 2020;Hendrikx et al., 2022).
The functional data was co-registered to the anatomical space using AFNI (afni.nimh.nih.gov;Cox, 1996) as previously described (Harvey et al., 2020;Hendrikx et al., 2022;Paul et al., 2022;van Ackooij et al., 2022).A single transformation matrix was constructed, incorporating all the steps from the raw data to the cortical surface model to reduce the number of interpolation steps to one.For the fMRI data, we first applied motion correction to the functional data (3dvolreg).For the timing scans, we also applied motion correction to the images that were acquired using opposing phase-encoding direction, then determined the distortion transformation between these and the functional runs (3dQwarp) to correct for spatial distortions in the functional scans (3dNwarpApply).Then we determined the transformation that co-registers this functional data to the T1 with the same resolution, position and orientation as the functional data (3dvolreg).We finally determined the transformation from this T1 image to a higher resolution (1 mm isotropic) whole-brain T1 image (3dUnifize, 3dAllineate).We applied the product of all these transformations for every functional volume to transform our functional data to the whole-brain T1 anatomy.We repeated this for each fMRI session to transform all their data to the same anatomical space.
The resulting data was imported into Vistasoft's mrVista framework (github.com/vistalab/vistasoft;Vistalab, 2021) for neural response modelling.For timing response data, we identified the parts of each E. Hendrikx et al. scanning run where each stimulus configuration was presented and averaged these fMRI responses together across all runs and sessions.For numerosity response data, we averaged all cycles within all runs together.We also separately averaged data from odd and even runs.

Neural response models
Our modeling approach follows a pRF modeling procedure (Wandell et al., 2007) as described previously (Harvey et al., 2020(Harvey et al., , 2015;;Harvey and Dumoulin, 2017a;Hendrikx et al., 2022;Paul et al., 2022).For each voxel, we find the parametric response function that most closely predicts the measured fMRI time course for that voxel.We repeat this for both the numerosity and timing response data sets.For numerosity response data, we analyzed the average response across all stimulus configurations with a single neural response model, as we have previously shown these produce very similar responses (Harvey et al., 2013;Harvey andDumoulin, 2017b, 2017a).For timing response data, we captured the responses to all stimulus configurations with a single model that predicts responses to changes in both event duration and period (Harvey et al., 2020;Hendrikx et al., 2022).
For a large set of candidate combinations of free parameters in each response model, we evaluated the response function at each presented numerosity or timing to predict the neural response time course this response function would produce when presented with our stimuli (Fig. 1, top rows).These neural response predictions were then convolved with a standard hemodynamic response function (HRF) to create a predicted fMRI response for each parameter combination (Fig. 1, bottom rows).These candidate fMRI response predictions are then compared to the recorded responses by computing their correlations (variance explained, R 2 ).
Having determined the parameters that best predict each voxel's responses, the resulting response function throughout the brain was then used to refit HRF parameters for each participant (Harvey and Dumoulin, 2011) because HRF parameters vary considerably between participants but relatively little between brain areas and recording sessions in a single participant (Handwerker et al., 2004).Having determined these participant-specific HRF parameters, the parameters of the response function were then refit following the same procedure.The response function parameters that most closely predict the fMRI time course of each voxel are then used for further analysis, together with the response variance explained by their predictions.

Numerosity response model
The best-fitting response function used to predict responses to numerosity is a logarithmic Gaussian function as described by Eq. ( 1) (Harvey et al., 2013;Harvey and Dumoulin, 2017a).
where Numerosity pref is the preferred numerosity, i.e., the mean of the logarithmic Gaussian response function, and σ is the range of numerosities this neuronal population is responsive to, i.e., the spread around the mean (standard deviation) in logarithmic numerosity space.These two free parameters allow an exhaustive high resolution search of their combinations.
After fitting this model, it is possible that some voxels are best described by a Gaussian function with a preferred numerosity outside of the presented stimulus range (i.e.below 1 or above 7).Similarly, some voxels' responses are poorly fit by the response model, with a variance explained less than or equal to 0.3, which is observable by chance with a probability above 0.05 (Harvey and Dumoulin, 2017a).In these cases it is not possible to accurately determine the voxel's preferred numerosity, or be sure that it shows a convincingly numerosity-tuned response.We therefore excluded such voxels from further analysis.Resulting parameters of this fitting procedure are given in Supplementary Fig. 2.

Timing response model
The best-fitting response function used to predict responses to event timing is a 2-dimensional anisotropic Gaussian function as described by Eqs. ( 2)-( 4) (Harvey et al., 2020;Hendrikx et al., 2022): where Duration pref and Period pref are the preferred duration and period, respectively, i.e., the mean of the Gaussian response function; σ maj and σ min are the range of timings this neuronal population is responsive to, i.
e., the spread around the mean (standard deviation); θ is the angulation of the major axis; and expFreq is the compressive exponent on frequency.This compressive exponent accounts for the sub-additive accumulation of fMRI responses with increasing event frequency.
Because of the large number of free parameters in this model, it is not possible to perform a high resolution exhaustive search over all parameter combinations.We therefore performed a low resolution exhaustive search followed by a gradient descent between the bestfitting parameter combination and its neighbors.
We predict a neural response with the amplitude described by Eqs. ( 2)-( 4) at the offset of every event.This is important because event frequency varies inversely with changes in event period in our stimuli, but responses to multiple events are aggregated due to fMRI's low temporal resolution.We tie this to the event's offset because the timing of the event is only known at its offset.
After fitting this model, it is possible that some voxels are best described by a Gaussian function with a preferred duration or period outside of the presented stimulus range (i.e.below 50 ms or above 1000 ms).Similarly, some voxels' responses are poorly fit by the response model, with a variance explained less than or equal to 0.2.Here we use a lower threshold because this 2-dimensional response function is sampled less thoroughly by our stimulus sequence and fit to a longer fMRI time series (including all 4 stimulus configurations).Furthermore, we have previously shown (Hendrikx et al., 2022) that voxels with a variance explained above 0.2 consistently allow us to distinguish between timing-tuned and monotonic responses, while those with a variance explained below 0.2 do not.Therefore, in such voxels it is not possible to accurately determine the voxel's preferred timing, or be sure that it shows a convincingly timing-tuned response.We therefore excluded such voxels from further analysis.Resulting parameters of this fitting procedure are given in Supplementary Fig. 3.

Regions of interest definitions
The borders of the numerosity and timing maps were defined based on the goodness of response model fits and the progression of response preferences across the cortical surface.Clusters of high model fits in cortical surface locations that were consistent across participants (Harvey et al., 2020;Harvey and Dumoulin, 2017a) were considered potential maps.These clusters could be split up into multiple neighboring maps if they contained multiple topographic progressions (e.g., timing map TTOP and TTOA or numerosity maps NPCM and NCPI, which are often touching).We named the numerosity and timing maps after their anatomical locations, following extrastriate visual field map naming conventions (Wandell et al., 2005).We preceded these locations with 'N' for numerosity maps and 'T' for timing maps (Harvey et al., 2020;Harvey and Dumoulin, 2017a) (example quantity maps: Fig. 2).Moving from posterior to anterior and inferior to superior, NLO (for 'numerosity temporo-occipital') lay in the lateral occipital sulcus.This small cluster was not characterized in our previous study (Harvey and Dumoulin, 2017a), as it was not clear identifiable in all participants in that study.We include it here because it is more consistently found in the current data set and lay close to a cluster of timing responsive voxels that we have previously described (Harvey et al., 2020) so is more important to our current research questions.NTO lay at the lateral temporo-occipital junction, between the inferior temporal and lateral occipital sulci, posterior-superior to the preoccipital notch.NPO lay superior to the parieto-occipital sulcus, medial to the posterior IPS.NPCI, NPCM, and NPCS (previously called NPC3, NPC2 and NPC1 in Harvey and Dumoulin, 2017a) lay in and around the postcentral sulcus.NPCI and NPCM lay in the postcentral sulcus, inferior and superior (respectively) to its junction with the IPS.NPCS lay immediately posterior to the superior postcentral sulcus and medial to the anterior IPS.NFI lay at the junction of the precentral and inferior frontal sulci: this numerosity map was not characterized in our previous study (Harvey and Dumoulin, 2017a).NFS (previously called NF in Harvey and Dumoulin, 2017a) lay at the junction of the precentral and superior frontal sulci.
Moving from posterior to anterior and inferior to superior, TLO ('timing lateral occipital') lay in the lateral occipital sulcus.TTOP and TTOA lay on the inferior lateral boundary of the temporal and occipital lobes, in the posterior inferior temporal sulcus.TPO lay superior to the parieto-occipital sulcus, medial to the posterior IPS.TLS lay in the posterior lateral sulcus.TPCI lay in the inferior postcentral sulcus.TPCM lay in the middle of the postcentral sulcus near its junction with the IPS.TPCS lay immediately posterior to the superior postcentral sulcus and medial to the anterior IPS.TFI lay at the junction of the precentral and inferior frontal sulci.TFS lay at the junction of the precentral and superior frontal sulci.

Overlap of numerosity and timing maps
The following new methods were used to quantify the overlap between numerosity and timing maps in this study.For each numerosity and timing map, we computed the proportion of responsive voxels that were also a responsive voxel in a map for the other quantity.As not all hemispheres have all maps, 119 numerosity map and 155 timing map examples were used.Additionally, we calculated the proportion of overlapping voxels that would be expected by chance: the proportion of all voxels within the hemisphere that are part of a map for the other quantity and convincingly responds to this other quantity.
First, for each map we assessed whether the observed proportion of overlapping voxels for the other quantity was greater than expected by chance.We performed one-sided approximate Wilcoxon signed-rank test for each map location, because the differences were not normally distributed in all maps (Shapiro-Wilk test).Here, we paired the observed proportion of overlap in the map example in each hemisphere with the chance overlap in the same hemisphere.False discovery rate correction (FDR correction Benjamini and Hochberg, 1995) was applied to the resulting probabilities as this test was repeated for all map locations, first over all numerosity maps and separately over all timing maps.Here, we calculated the effect size as r = Z / ̅̅ ̅ n √ , where n is the number of pairs.
Second, to assess whether the proportion of overlap differs between the maps, the proportion of overlap was compared using a 3-factor ANOVA (factors: participant, map location, and hemisphere) using MATLAB's anovan function.We include hemisphere as a factor because we have previously found differences between hemisphere in numerosity (Harvey and Dumoulin, 2017a) and timing (Harvey et al., 2020) response function parameters.Post-hoc multiple comparisons were then performed using MATLAB's multcompare function with Tukey-Kramer correction for multiple comparisons, because the maps differed significantly and the residuals of this ANOVA were normally distributed (Shapiro-Wilk).We repeated these analyses for any voxels responding to the other quantity that did not have to be part of a specific map for the other quantity, or indeed any region of interest: this avoided dependance on our defined regions of interest when characterizing the overlap of responsive areas (Supplementary Fig. 4; N numerosity maps = 119; N timing maps = 155).Furthermore, we also repeated these analyses limiting the overlapping voxels for the other quantity to the map of the other quantity that most overlaps each timing and numerosity map and found similar results (Supplementary Fig. 5; N numerosity maps = 119; N timing maps = 155).
Third, we tested whether the proportion of overlap changed over the posterior-anterior and inferior-superior axes of the brain.Here, we first transformed each participant's anatomical MRI data, together with each map example's center, to the Talairach N27 surface using AFNI's 3dAllineate and 3dNwarpApply tools.We then used an ANCOVA to test whether the proportion of overlap (in the participant's native space) was predicted by the y and z Talairach coordinates of the map examples (continuous factors), taking participant and hemisphere into account as additional categorical factors.We applied FDR correction (Benjamini and Hochberg, 1995) to the resulting probabilities as this test was repeated for numerosity maps and timing maps.
Fourth, we assessed whether the proportion of overlapping voxels had a relation to the total amount of voxels responding to the quantity for each map example.We used an ANCOVA to test whether the proportion of overlap was predicted by map size (continuous factor), taking participant and hemisphere into account as additional categorical factors.We applied FDR correction (Benjamini and Hochberg, 1995) to the resulting probabilities as this test was repeated for numerosity maps and timing maps.
Fifth, we assessed whether the proportion of overlap between numerosity and timing maps within each map location differed from the proportion of overlap during repeated measures of the same quantity (Supplementary Fig. 6).For repeated measures we split the data into even and odd scans.For each numerosity and timing map example, we computed the proportion of responsive voxels that were also a responsive voxel in a map for the same quantity in the other half of the data.We took the average proportion of the two halves for further analyses.As this proportion of overlap is only based on splits of the data, the measurements are more noisy.To account for this discrepancy in noise, we compared these repeated measures of the proportions to the proportion of overlap of numerosity and timing maps also using split halves.Here, we performed the analysis described above, but quantifying the proportion of overlap between (1) numerosity maps in the odd scans and timing maps in the even scans; and (2) numerosity maps in the even scans and timing maps in the odd scans.We then averaged these two proportions.We compared the proportion of overlap within each map location using a two-sided Wilcoxon signed-rank test, because the difference was not normally distributed (Shapiro-Wilk) (N numerosity = 119; N timing maps = 156).FDR correction (Benjamini and Hochberg, 1995) was applied to the probability of observed differences, first over all numerosity maps and separately over all timing maps.Here, we calculated the effect size as r = Z / ̅̅ ̅ n √ , where n is the number of pairs.

Correlation
We also tested whether numerosity and timing preferences within the overlapping areas were related.For each map we looked at voxels convincingly responding to the other quantity, that were a part of a map for the other quantity.There are several complexities to this question.First, Gaussian tuning functions for numerosity and timing are in logarithmic and linear spaces, respectively.We therefore use Spearman's rank (non-parametric) correlations.Second, if two completely unrelated topographic progressions overlap, a correlation would nevertheless be expected (whether negative or positive) unless the progressions are perfectly orthogonal.However, if the progressions are unrelated the sign of this correlation should not be consistent across map examples.We E. Hendrikx et al. therefore assess the distribution of correlation coefficients across map examples rather than the significance of the correlation within individual map examples.We limit this to map examples with at least 10 overlapping voxels, as correlation coefficients across very small numbers of measurements are not meaningful.As TLS has no analogous numerosity map and only shows overlap for one map example, we did not consider this a reliable correlation distribution.Therefore we excluded this datapoint from all following analyses.As a result, 97 numerosity map and 107 timing map examples were used here.
For each numerosity map and timing map example we found all overlapping maps for the other quantity and their responsive voxels.We then determined the Spearman's rank correlation coefficient between numerosity preferences and both the duration and period preferences within that overlapping voxel set.For each map location (e.g.all examples of NPCS) we took the correlation coefficients from all map examples.We found that these were not normally distributed in each map location (Shapiro-Wilk test).We therefore performed two-sided onesample Wilcoxon signed-rank tests to assess whether the correlation between quantity preferences within each of the map locations was significantly different from 0 for all map locations with more than one map example.We FDR corrected (Benjamini and Hochberg, 1995) the resulting probabilities across all map locations and across comparisons between numerosity preferences and both duration and period preferences.To maximize this test's sensitivity to any consistent relationship, we also repeated this analysis using the correlation coefficients from all map examples in all locations together.Here, we calculated the effect size as r = Z / ̅̅ ̅ n √ , where n is the number of included map examples.Second, to assess whether the correlation between preferred values differs between the maps, the correlation coefficients were compared using a 3-factor ANOVA (factors: participant, map location, and hemisphere) using MATLAB's anovan function.Post-hoc multiple comparisons were then performed using Dunn's test with Holm-Šidák correction for multiple comparisons (Sidak, 1967), because the maps differed significantly and the residuals of this ANOVA were not normally distributed (Shapiro-Wilk).
Note that it is possible that a single map for one quantity would overlap with multiple maps of another quantity.This could show a strong correlation to one of these maps, but this could be hidden by correlations with the opposite sign to other maps.Therefore, we also performed these analyses focusing on each map's relation with its most overlapping map for the other quantity (N numerosity maps = 97 and N timing maps = 107; Supplementary Fig. 7).
It is possible that measures of preferred quantity are less precise in map locations with a smaller range of this quantity.This would cause correlations in maps with a smaller range of quantities to be lower.To test whether the range of the quantities had a (non-parametric) relationship to the strength of the correlation, we used partial Spearman correlations.We assessed the relation between absolute correlation between the preferred values of the map examples with the interquartile range of log(preferred numerosity) and preferred timing (Supplementary Results Text 1).We applied FDR correction (Benjamini and Hochberg, 1995) to the resulting probabilities as this test was repeated for duration and period preferences in numerosity maps and timing maps.
The distribution of correlation coefficients across map examples might not differ significantly from zero if data quality for either response model was too poor to show repeatable quantity preferences across repeated measurements.To exclude this possibility, we performed the same analysis comparing preferred numerosity estimates from odd and even scan runs, comparing preferred duration estimates from odd and even scan runs, and comparing preferred period estimates from odd and even scan runs.Notably, each of these analyses uses half as much data as in the main analysis, so should be less sensitive to correlation coefficients differing significantly from zero.We use the same method for selecting voxels as described above: we found all overlapping maps for the same quantity and their responsive voxels in the other half of the data.
We then compared the consistency of correlations between preferences during repeated measures with the consistency of correlations between numerosity and timing preferences.To ensure a fair comparison, we used numerosity and timing preference estimates from odd and even runs in both cases.We computed the correlations between the odd numbered numerosity and even numbered timing runs and averaged this with the correlations between even numbered numerosity and odd numbered timing runs, using Z-transformations.If there was only data available for one of these comparison, this comparison was taken as the average (N numerosity maps = 93; N timing maps = 102; Supplementary Fig. 8).Note that NLO and TLO do not have more than one map example with overlapping numerosity and timing preferences and are thus excluded from these analyses (like TLS before).We compared the correlations between numerosity and timing preferences, with the correlations between repeated measures for each map location using a two-sided Wilcoxon signed-rank test, because the difference was not normally distributed (Shapiro-Wilk).FDR correction (Benjamini and Hochberg, 1995) was applied to the probability of observed differences, first over all numerosity maps and separately over all timing maps.Here, we calculated the effect size as r = Z / ̅̅ ̅ n √ , where n is the number of pairs.Repeated measures of maps have more overlapping voxels than numerosity and timing maps in halves of the data.The previous analysis may therefore bias the repeated measures to have a higher data quality and consequently higher, less noisy, correlations.Therefore, the above analysis was also repeated selecting only voxels that were included in all compared maps.This comparison has less datapoints, as not all the compared maps have at least 10 overlapping voxels (N numerosity maps = 89; N timing maps = 96; Supplementary Fig. 8).

Topographic progressions
A significantly non-zero distribution of correlations between numerosity and timing preferences could be absent because the topographic progressions of the preferred numerosities and timings within overlapping maps were not closely aligned, but still had a consistent relationship between quantity preference progressions (e.g., if their progressions are orthogonal to each other).This would still be a meaningful relationship, but could produce a distribution of correlation coefficients that is not significantly different from zero.Alternatively, the topographic progressions could be in completely unrelated directions in different map examples, so have no relationship at all.To distinguish these possibilities, we investigated the relationships between the direction of topographic progressions of numerosity and timing preferences.
For each voxel in the overlapping set, we determined the direction of greatest change in numerosity and timing preferences.This direction also considered changes among responsive voxels outside the overlapping region.Each map's topography was determined using a method based on Guo and colleagues (Guo et al., 2012).Gradient vectors g ij of the preferred values were calculated between each of the voxels within a map.This was done by firstly calculating the difference of the preferred quantity (log (numerosity), duration or period) between every source voxel i and every target voxel j (v i , v j ) within the map.This difference was normalized by cortical surface distance (denoted by subscript) between the voxels to take into account that larger cortical distances should have larger changes in quantity preferences (Eq.( 5)).
Next, to give the gradient vector a direction running along the cortical surface from the source voxel (i.e.orthogonal to the surface normal), accounting for the folded surface of the human brain, we adapted Guo and colleagues' method.We first determined the vertices falling on the shortest path (X) between the source (v i ) and target voxel (v j ).We reasoned that longer distances along this path would give a E. Hendrikx et al. more precise estimate of its direction, as long as the path did not cross the surface normal (i.e. the cortical surface did not fold back upon itself).We therefore stepped along the path calculating the cross product of the surface normal of source voxel v i (n i ) and the straight-line direction from the source voxel v i to the current voxel along the path (v x ) (i.e. the direction travelled along the cortical surface to follow this path).We compared the surface direction of this voxel (v x ) with the surface direction of the first step along path X (i.e., the step from v i to v 1 ) (Eq. ( 6)).
where x ranges from 1 until the length of the shortest path between v i and v j (X); the vector v i to v 1 is the initial step along path X; and n i is the surface normal at voxel v i .Here, × is used in the vector context, so denotes a cross product, and ∠ is the angle between the two terms in brackets.
We selected voxels (v y ) for which the surface direction differed less than 90 • from the surface direction of the first step (suggesting the path had not folded back on itself) (Eq.( 7)).
We took the final straight-line direction of the voxel furthest along this path for which this was true as the direction or the vector from the source to the target voxel, and calculated its unit vector (Eq.( 8)).
where v max is max(v y ) We multiplied this direction with the distance we had already determined to give the vector from source voxel v i to target voxel v j (Eq.( 9)).
Here, × is used in the scalar context so denotes a scalar multiplication (although the result is a vector).
We repeated this for all target voxels then averaged the resulting vectors to give the direction of greatest increase in quantity (Eq.( 10)).
Finally, we took the cross product of the resulting vector and the surface normal to give the direction along the cortical surface at the source voxel.
For each voxel v i we determined the angle between these vectors for numerosity and duration preferences, and numerosity and period preferences.Regardless of whether we were looking in a numerosity or timing map, we determined the angle of the timing progression relative to the numerosity progression, such that the angle of the numerosity progression was zero.In the left hemisphere maps we reversed the sign of the angle to account for the reversal of handedness between the hemispheres (i.e. the hemispheres are mirror reversed).Therefore, all angles are given relative to medial-lateral rather than the left-right directions, though this choice does not change the directional effects we find (Supplementary Fig. 9; Supplementary Table 22).
We then computed the circular mean of all angles within a overlapping set of voxels using circ_mean from the circular statistic toolbox for MATLAB (Berens, 2009).Note that not all maps show overlap with other maps.As not all map examples have any overlapping map for the other quantity or have fewer than 10 overlapping voxels, 97 numerosity map and 107 timing map examples were used in this analysis.
To assess whether the angles had a consistent distribution within the maps, we applied a Bayesian analysis of circular uniformity using the BayesCircIsotropy package in R (Mulder and Klugkist, 2021).Here we test the null hypothesis of uniformly distributed angles against an alternative hypothesis of a von Mises distribution with a conjugate prior of p(κ) ∝I 0 (κ) − 1 , where κ ∈ R + is a concentration parameter and I 0 (⋅) is the modified Bessel function of the first kind and order zero.This is a recommended prior for a low expected concentration within this distribution (Mulder and Klugkist, 2021).This is tested against the null hypothesis is a uniform circular distribution of vector angles.
Again here, it is possible that a single map for one quantity would overlap with multiple maps of another quantity, so a strong topographic alignment in one map could be hidden different relationships in other maps.Therefore, we again repeated this analysis using each map's relation with only its most overlapping map for the other quantity (N numerosity maps = 97 and N timing maps = 107; Supplementary Fig. 10).
Then, we aimed to quantify potential increases or decreases in consistency of the topographic alignments over the posterior-anterior and inferior-superior axes of the brain.To this end, we computed the partial Spearman correlations between the Bayes factors; the mean y and mean z Talairach coordinates (Supplementary Results Text 2).We applied FDR correction (Benjamini and Hochberg, 1995) to the resulting probabilities as this test was repeated for numerosity maps and timing maps.We also tested the possibility that significant clustering of vector angles could be found only in maps with larger overlapping areas and so potentially less noisy angle measures.To look for this effect across individual map examples, we tested whether the number of overlapping voxels within each map example significantly predicted the standard deviation of the vector angles within that map example (Supplementary Results Text 3).Here we used an ANCOVA with number of overlapping voxels, participant, and hemisphere as factors.FDR correction (Benjamini and Hochberg, 1995) was applied to the probability of observed differences as this test was repeated for duration and period preferences in numerosity maps and timing maps.To look for a similar effect across map locations, we tested whether the map location's average number of overlapping voxels was significantly (Spearman) correlated with the Bayes factor of the von Mises distribution associated with this map location (Supplementary Results Text 3).FDR correction (Benjamini and Hochberg, 1995) was applied to the probability of observed differences as this test was repeated for duration and period preferences in numerosity maps and timing maps.
We then tested the possibility that significant clustering of vector angles could be found only in maps with a larger range of preferred quantities in the overlapping areas, potentially allowing less noisy angle measures.Therefore, we tested whether the map location's median interquartile ranges for the preferred quantities in the area of overlap were significantly (Spearman) partially correlated with the Bayes factor of the von Mises distribution associated with this map location (Supplementary Results Text 4).FDR correction (Benjamini and Hochberg, 1995) was applied to the probability of observed differences as this test was repeated for duration and period preferences in numerosity maps and timing maps.
As with the correlation analyses, the distribution of vector difference directions across map examples may not differ from uniform if data quality for either response model was too poor to show repeatable directions of progression across repeated measurements.To exclude this possibility, we performed the same analysis comparing preferred numerosity estimates from odd and even scan runs, comparing preferred duration estimates from odd and even scan runs, and comparing preferred period estimates from odd and even scan runs.We use the same method for selecting voxels described above: we found all overlapping maps for the same quantity and their responsive voxels in the other half of the data.We computed the average angles between the odd and even runs using circ_mean.If there was only data available for one of the runs, this run was taken as the average (N numerosity maps = 117; N timing maps = 148; Supplementary Fig. 11).Notably, each of these analyses uses half as much data as in the main analysis, so should be less sensitive to non-uniform distributions of vector angles than the main analysis.
We then compared the consistency of the angles between topographic progressions found in repeated measures of the same quantity E. Hendrikx et al. preferences to that found between numerosity and timing preferences.To ensure a fair comparison, we again used the numerosity and timing preference estimates from odd and even runs in both cases.We computed one angle between the progressions estimated from odd numbered numerosity runs and even numbered timing runs and averaged this with the angle between progressions estimated from even numbered numerosity runs and odd numbered timing runs, using circ_mean.If there was only data available for one of these comparisons, this comparison was taken as the average (N numerosity maps = 95; N timing maps = 102; Supplementary Fig. 12).We compared the Bayes factors supporting a von Mises distribution rather than a uniform distribution of angles between numerosity and timing preferences, with those of repeated measures.Here we used a two-sided Wilcoxon signed-rank test, because the difference between these Bayes factors was not normally distributed (Shapiro-Wilk).Note that this comparison is only possible using all map locations together (rather than for each map location separately) as each quantity map location only yields one Bayes factor (Supplementary Table 22).FDR correction (Benjamini and Hochberg, 1995) was applied to the probability of observed differences, first over all numerosity maps and separately over all timing maps.Here, we calculated the effect size as r = Z/ ̅̅ ̅ n √ , where n is the number of pairs.Repeated measures of maps for a single quantity have more partially overlapping maps and more overlapping voxels within these maps than comparisons between numerosity and timing maps.The above analysis may be affected by the resulting differences in data quality and consequently the precision of angles estimates.To ensure that this analysis wasn't therefore biased in favor of repeated measures, we also repeated it selecting only voxels that were included in all compared maps (so the same voxel selection in all comparisons).This comparison has less datapoints, as not all the compared maps have at least 10 overlapping voxels (N numerosity maps = 89, N timing maps = 96; Supplementary Figs.13-14; Supplementary Table 22).

Increasing overlap between numerosity and timing maps along inferior-superior hierarchies
If numerosity and timing responses were transformed into common magnitude responses, their responses should overlap.We looked at whether the overlap between the numerosity and timing maps was significantly greater than expected by chance (i.e. if their positions were unrelated, see Methods).We found that for most maps the proportion of overlapping voxels with a map for the other quantity was greater than chance (p < 0.01; Fig. 3; Table 1).Exceptions were the most posterior numerosity (NLO) and timing (TLO) maps; an inferior timing map that does not have a numerosity map in a similar location (TLS); and timing map TTOA.Numerosity map NTO is posterior to TTOA, largely overlapping with the adjacent timing map TTOP.
Next, we tested whether the proportion of overlap between the numerosity and timing maps changed along their posterior-anterior and inferior-superior hierarchies (Harvey et al., 2020;Harvey and Dumoulin, 2017a;Hendrikx et al., 2022).We used two ANOVAs (one for numerosity and one for timing maps) with map location, hemisphere and participant as factors.
We found similar results when repeating these analyses without limiting the overlap to voxels within a map for the other quantity Fig. 3. Hierarchical increase in the proportion of overlap between numerosity and timing maps.(A) Moving from posterior and inferior (left of the panel) to anterior and superior (right of the panel) numerosity maps, the proportion of voxels within a numerosity map that also show a timing-tuned response increases.Except for the most posterior numerosity map (NLO), the proportion of overlap (colored shapes) is significantly greater (stars below the axis) than expected by chance (colored lines).(B) Moving from posterior and inferior (left) to anterior and superior (right) timing maps, the proportion of voxels within a timing map that also show a numerosity-tuned response increases.Except for the most posterior and most inferior timing maps (TLO, TTOA, TLS), the proportion of overlap is significantly greater than expected by chance.Different colors show results from different participants, circles show left hemisphere maps and diamonds show right hemisphere maps.Squares show the mean proportion of overlap across all examples of the same map, error bars show the standard error.Brackets and stars show significant differences in multiple comparisons with maps grouped across hemisphere: all brackets to the left of the star are significantly different from all brackets to the right of the star (p < 0.05).(For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)E. Hendrikx et al. (Supplementary Fig. 4; Supplementary Tables 3-5) or while limiting the overlap to the most overlapping map for the other quantity (Supplementary Fig. 5; Supplementary Tables 6-7).
We have previously shown that anterior and superior numerosity (Harvey and Dumoulin, 2017a) and timing maps (Harvey et al., 2020) tend to be larger than posterior and inferior maps.If two nearby large maps have a higher chance of overlapping than two nearby small maps, could the observed increase in overlap reflect this increase in map size?Using ANCOVAs taking participant (numerosity: F (7109) = 3.07, p = 0.005, ƞ p 2 = 0.16; timing: While these effect sizes are smaller than those for map z coordinate, map overlap differences may partially reflect map size differences.
If numerosity and timing were transformed into a common representation, the overlap between numerosity and timing maps should be similar to the overlap between repeated measures of numerosity maps or of timing maps.However, the proportion of overlap between maps for different quantities was significantly less than the proportion of overlap between repeated measures of all numerosity maps, and repeated measures of all timing maps except for TPO, TPCM and TPCS (which show no significant difference) (Supplementary Fig. 6; Supplementary Table 8).Therefore, numerosity and timing maps are different structures despite partially overlapping.

Correlations between numerosity and timing preferences only in maps around intraparietal sulcus (IPS)
If numerosity and timing responses were transformed into common magnitude responses, numerosity and timing preferences within the maps should be correlated consistently and as strongly as repeated measures of either quantity preference.We first looked at the distribution of correlation coefficients between numerosity and timing preferences.The correlation within each overlapping area reflects data quality as well as any relationship between response preferences.Therefore, if numerosity and timing preferences are related in a particular map location, their correlation might reach significance in only some examples but its direction (sign) should be consistent across examples of the same map.Conversely, if there is no relationship, the correlation should not have a consistent direction, although it might reach significance in some examples of a particular map location.We therefore focused on the repeatability of the sign of correlations across participants: Is the distribution of correlation coefficients (measured from all examples of the same map location) significantly different from zero?
In maps surrounding the IPS (NPO, NPCM, TPO), the distribution of correlation coefficients was significantly larger than zero (Fig 4; Table 2).The timing map TPCM (near numerosity map NPCM) showed a distribution of correlation coefficients that just failed to reach significance after correction for multiple comparisons.This difference did not reach significance in other maps, although numerosity map NFS and timing map TFS show trends towards negative correlations between numerosity and duration preferences.We also tested the distribution of correlation coefficients from all map examples to test for a global relationship between these quantity processing networks, and found no significant difference from zero.Therefore, relationships between numerosity and timing preferences appear to be limited to some intraparietal areas only.
Limiting correlations to the most overlapping map only (to exclude possible effects of different correlations in different maps) produced similar results (Supplementary Fig. 7; Supplementary Tables 13-17).While the ranges of numerosity and timing preferences differ between map locations (Harvey et al., 2020;Harvey and Dumoulin, 2017a), the (absolute) correlation coefficients showed no significant relationship to the range of quantity preferences in compared maps (see Inline Supplementary Results Text 1).
To demonstrate the data quality in these map locations was sufficient to detect correlations between quantity preferences, we repeated these analyses on quantity preferences estimated from repeated measures of numerosity and timing mapping data (i.e.split halves of the data, for the same quantity).The distribution of correlation coefficients was significantly above zero in all but one map location (TFI for repeated measures of duration preferences) (Supplementary Fig. 8; Supplementary  Table 18).Therefore, although estimated from only half the data, repeated measures of preferences for the same quantity were more consistently related than measures of different quantities.Furthermore, the distribution of correlations from repeated measures of the same quantity was significantly higher than the correlations between preferred numerosity and timing except for numerosity and duration preferences in TPO, TTOA and TFI; and numerosity and period preferences in TPO (where there was no significant difference) (Supplementary Fig. 8; Supplementary Table 19).This generally contradicts the proposal that numerosity and timing are transformed onto a common representation, even in the limited overlapping area.However, as for the proportion of overlap, TPO shows similarly close relationships between different quantity preferences and repeated measures of the same quantity, supporting a common representation in this specific area.
Because numerosity and timing maps are not the same structures, they have fewer overlapping voxels than repeated measures of the same quantity map, giving less statistical power.We therefore repeated the previous comparisons using only voxels within both the repeated measure overlaps and the different quantity overlap: the same set for both correlations (Supplementary Fig. 8; Supplementary Tables 20 -21).This smaller voxel selection produced similar results, although fewer differences reached significance.

Consistent angles between topographic progressions in some numerosity and timing maps
This lack of correlation between different quantity preferences in overlapping areas may reflect 1) noisy preferred quantity estimates (unlikely given the consistent correlations between repeated measures); 2) a lack of consistent relationships between preferences; or 3) a consistent angle between quantity preference progressions that does not result in a strong correlation, e.g., the progressions are orthogonal to each other.
Therefore, we asked whether the direction of topographic progression of numerosity and timing preferences was consistently related across examples of each numerosity and timing map (see Methods and Fig. 5 for a schematic overview).We computed the angle between these progressions in a circular space, where 0 • indicates a perfect alignment of numerosity and timing progressions; 180 • or − 180 • indicates progressions in opposing directions; 90 • and − 90 • indicate increases in the two possible orthogonal directions.
For some quantity maps, Bayesian analyses revealed that angles between numerosity and timing progressions over map examples supported von Mises (circular Gaussian) distributions more than uniform distributions, i.e. angle between progressions across examples somewhat consistent between examples of the same map location.In TPCS and NPCS the evidence for consistent angles between the numerosity and period progressions was even considered very (BF 30-100) or extremely strong (BF>100), respectively.Overall, the consistent relationships did not follow a clear hierarchy (see Supplementary Results Text 2).
The mean angle between quantity progressions centered around 0 • for some more posterior maps (NTO, NPO, TTOP, TPO.Fig. 6), but around ±180 • for superior parietal and frontal maps (NPCS, NFS, TPCS, TFS.Fig. 6A and B).Similar results were found when using left-right rather than medial-lateral directions (Supplementary Fig. 9) and when considering only the most overlapping map rather than all overlapping maps (Supplementary Fig. 10).The standard deviations and Bayes factors in each map location were not significantly predicted by the number of overlapping voxels or the range of quantity preferences in each map location (see Supplementary Results Texts 3-4).
To demonstrate that data quality was sufficient to detect consistent angles between topographic progressions, we repeated this analysis between repeated measures (i.e., split halves) of numerosity and timing preferences (Supplementary Fig. 11).These preferences are noisier, being derived from only half the data.This analysis supported von Mises distributions consistently centered near 0 • in all map locations except for preferred duration estimates in TPCI and TFI.Therefore, where angles between progressions were inconsistent, this reflected variable relationships between overlapping maps rather than insufficient data quality.Furthermore, the consistency of angles between progressions from repeated measures of the same quantity was significantly higher than the consistency of angles between progressions of numerosity and timing preferences (Supplementary Fig. 12; Supplementary Table 22).Note that, since we only have one measure of consistency per map location, we cannot make this comparison separately per map location.
Because numerosity and timing maps are not the same structures, they have fewer and smaller areas of overlap than repeated measures of the same quantity map, potentially giving less statistical power.We therefore repeated the previous comparisons using only voxels within both the repeated measure overlaps and the different quantity overlap: the same set for both angles between topographies.This smaller voxel selection produced similar results, although fewer differences reached significance (Supplementary Figs. 13 -14; Supplementary Table 22).Overall, these results indicate that the topographic progressions of numerosity and timing preferences are not consistently related throughout the brain, though they may show some relationships (both negative and positive alignments) in some maps.This generally contradicts the proposal that numerosity and timing are transformed onto a common representation.

Discussion
Behavioral and neural commonalities and interactions between quantities have led to the influential hypothesis that the representations  of different quantities may converge into a single, shared representation of magnitude (Walsh, 2003).We therefore examined the overlap between cortical areas showing numerosity and timing selective responses; the relationship between numerosity and timing preferences (and their topographic progressions) within these overlapping areas; and how these measures changed over their respective cortical processing hierarchies.We find that the earliest (most posterior) numerosity and timing maps showed little to no overlap.After this, areas responding to numerosity and timing overlap more than expected by chance and are increasingly brought together, mostly over their inferior-superior hierarchies.Within the overlapping regions, we found significant correlations between numerosity and timing preferences only in specific maps around the IPS.To investigate similarities in quantity preference progressions between numerosity and timing maps, we looked at the angle between these topographic progressions within the overlapping areas.In line with our correlation analyses, in some maps these quantity preference progressions were consistently aligned, but in some anterior superior areas the quantity preferences consistently progress in opposite directions.Overall, these results, suggest that neural responses to numerosity and timing are initially derived separately (Hendrikx et al., 2022;Paul et al., 2022).While these responses overlap later in the hierarchy, the lack of increasingly consistently correlated numerosity and timing preferences, together with a lack of increasingly aligned topographic progressions, is not consistent with a hierarchical transformation into a generalized and common quantity representation.
In previous neuroimaging studies (Dormal et al., 2012;Pinel et al., 2004;Skagerlund et al., 2016), common representations of different quantities have been proposed because of overlapping activations during quantity comparison tasks.Having established hierarchical processing of numerosity and timing responses (Harvey et al., 2020;Harvey and Dumoulin, 2017a), we can build on these earlier findings to show that the overlap between maps for different quantities increases over these quantity processing hierarchies.
Mapping quantity preferences further allows us to investigate the relationships between preferences within these overlapping areas.Transformation onto common quantity representations would predict increasingly related responses together with the increasing overlap.However, we do not see an increasingly close relationship between Fig. 6.Consistent angles between preferred numerosity and timing progressions in some maps (A) Angles between preferred numerosity and preferred duration progressions for each numerosity map, ordered from posterior (left) and inferior (bottom) to anterior (right) and superior (top).The angles support a von Mises distribution over a uniform distribution in some numerosity maps (NTO, NPO, NPCS, NFS) (Bayes factors BF 10 below map name).(B) As in (A), but for angle between preferred numerosity and preferred period progressions and supporting a von Mises distribution in NTO and NPCS.(C) Angles between preferred numerosity and preferred duration progressions for each timing map, supporting a von Mises distribution in the superior timing maps (TPO, TPCS, TFS).(D) As in (C), but for angle between preferred numerosity and preferred period progressions and supporting a von Mises distribution in TTOP and TPCS.Different colors show angles between quantity progressions from different participants, circles show left hemisphere maps and diamonds show right hemisphere maps.Squares show the circular mean of the angles across all examples of the same map location, error bars show the circular standard deviation.(For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)numerosity and timing preferences, which are correlated around the IPS: in a posterior parietal area (where numerosity map NPO and timing map TPO overlap), and in a postcentral overlapping area (where numerosity map NPCM overlaps with timing maps).Other areas did not show consistent correlations.Correlations between quantity preferences were generally lower and less consistent than correlations between repeated measures of preferences for one quantity.This suggests that these different quantities are not transformed into a common response.Indeed, when all maps were considered together, we observed no consistent correlations between numerosity and timing preferences because their correlations differed in strength and sign between maps.
Similar to relationships between numerosity and timing preferences we describe here, we have previously found consistent correlations between numerosity and object size preferences in the IPS (Harvey et al., 2015), though in that case we only looked in the IPS.This is striking because IPS has been speculated to be a potential location for magnitude-general representations (Bueti and Walsh, 2009).Although relations within the overlap of one small and specific area are indistinguishable from a common response to quantities, it is hard to draw strong conclusions about relationships.On the one hand, many nearby recording sites (those in NPO but outside TPO) respond to only numerosity.Furthermore, we expected a hierarchical emergence of a common response, in line with hierarchical changes in other properties of numerosity and timing responses (Harvey et al., 2020;Harvey and Dumoulin, 2017a;Hendrikx et al., 2022;Paul et al., 2022).The absence of similar relationships throughout the network does not support this assumption.On the other hand, the assumption of a linear hierarchy may be too simplistic.Potentially, parallel pathways each integrate quantity representations and one of these pathways indeed may lead to a common representation specifically in the overlapping area of NPO and TPO.
Previous results suggesting common quantity representations around IPS are complicated by the use of similar comparison tasks for both quantities (Dormal et al., 2012;Pinel et al., 2004;Skagerlund et al., 2016).These studies contrast responses during quantity comparisons to responses without these comparisons.They attribute common responses during both comparisons to a common representation of both quantities, rather than a common task being used for both quantities (DeWind et al., 2017;Harvey et al., 2017).Notably, the IPS is implicated in comparisons of quantities, symbolic numbers and non-quantitative features alike (Cohen Kadosh et al., 2005;Dormal et al., 2012;Fias et al., 2002;Göbel et al., 2004;Harvey et al., 2017;Pinel et al., 1999Pinel et al., , 2001Pinel et al., , 2004;;Shuman and Kanwisher, 2004), suggesting this may be a response to comparisons and not to quantities.
Here we avoided this issue by mapping neural quantity preferences in the absence of quantity comparisons.While we use a common color change detection task while mapping responses to both quantities, this task was ongoing throughout the whole experiment.This design distributes neural responses to our task evenly across all quantity magnitudes, while our analyses only looks at signals associated with specific quantity magnitudes.Therefore quantity preference estimates and model fits should be unaffected by our tasks.
Despite this, we still found large overlaps between responses to these different quantities, particularly in parietal and frontal areas implicated in comparisons and action planning (Sawamura et al., 2002).We therefore speculate that distinct sensory responses to different quantities may be brought together in these common areas to allow sharing of general comparison computations (DeWind et al., 2017) but largely without transforming these quantities into common representations.
Mapping quantity preferences reveals topographic quantity preference progressions.Assessing the relationships between these topographic progressions lets us investigate relationships between numerosity and timing preferences at a larger scale than voxel-to-voxel correlations.Transformation onto common quantity representations would predict increasingly aligned topographies as the overlap between the maps increases.However, we only find positive alignments for temporal occipital (NTO, TTOP) and parietal occipital (NPO, TPO) areas.These alignments to numerosity topography do not always hold for both duration and period.
Notably, some of the angulations between topographies we see in the superior parietal and frontal areas (NPCS, TPCS, NFS, TFS) have nonuniform distributions centered around 180 • : numerosity and timing preferences change in approximately opposite directions.This may be because we express timing responses in terms of period rather than frequency (1/period), so when period increases frequency decreases.Responses in more anterior timing maps increasingly depend on period or frequency rather than duration (Harvey et al., 2020).This opposite alignment therefore associates high numerosities with high frequencies rather than long durations.Different directions of alignments in different areas could potentially explain the different directions of the behavioral numerosity and timing interactions found in different studies (Javadi and Aichelburg, 2012;Lambrechts et al., 2013;Martin et al., 2017;Tsouli et al., 2019).
While the correlation and topographic alignment analyses may appear to test for similar relationships between numerosity and timing preferences, they sometimes give surprisingly different results.Not all areas showing consistent correlations also show consistent alignments or vice versa.Indeed, individual map examples with low but positive correlations between quantity preferences could show opposite directions of quantity progressions.These discrepancies are straightforwardly understood by considering that the correlation analyses only consider the overlapping areas.While the angles between topographic progressions can only be calculated within overlapping areas, these underlying directions of progression are based on entire quantity maps.
Even where different quantities are not clearly transformed into common responses, the observed overlapping responses to different quantities may be enough for neural populations responding to different quantities to interact and elicit behavioral interactions.Indeed, numerosity is represented in a distributed population code: information about numerosity 3 is not only held by neural populations with a preferred numerosity of 3, but also by the changing responses of neural populations with other numerosity preferences (Eger et al., 2009;Tudusciuc and Nieder, 2007).Therefore, all of the informative neural populations could interact in overlapping maps, rather than just the most responsive neural populations.A one-to-one correspondence between two quantities' preferences does not seem necessary to allow interactions.
Furthermore, we might expect that many maps with one-to-one correspondences between quantity preferences would lead to a strong, symmetric interaction of the quantities' perceptions, or even confusion between one quantity and one another.The lack of close correspondence that we observe here is consistent with the weaker, asymmetric behavioral interactions are often reported (Agrillo et al., 2010;Alards-Tomalin et al., 2016;Dormal et al., 2006;Lambrechts et al., 2013;Martin et al., 2017;Tsouli et al., 2019).
As we mapped numerosity and timing preferences in separate experiments, we did not investigate how one quantity affected responses to the other.Recent work shows that when numerosity and timing change simultaneously, their preferences change as a result of quantity interaction (Fortunato et al., 2023).Likewise, quantity preferences are affected by the context of recent states of the same quantity (Kristensen et al., 2021;Tsouli et al., 2022), even though the previously-viewed quantity activates a distinct neural population from the current quantity: adaptation to numerosity affects the whole numerosity map's response to all quantities.So any effect on quantity-selective responses may be sufficient to affect perception of that quantity.
Finally, it is possible that neural interactions occur before tuned responses emerge.Early visual areas monotonically increase their response amplitudes at the retinotopic location of the stimulus with numerosity, event duration and event frequency.It is unlikely that these similar early visual encoding schemes directly lead to a common representation of these quantities because the first areas showing tuned E. Hendrikx et al. responses to numerosity and timing are distinct (Hendrikx et al., 2022;Paul et al., 2022).Nevertheless, common early visual encoding schemes may underlie perceptual commonalities between these quantities, while perceptual interactions may result from interactions between simultaneous monotonic responses when both quantities vary in the same task.

Conclusion
Therefore, our results do not strongly support the most literal interpretation of A Theory Of Magnitude: the transformation of different quantities onto a common metric throughout cognitive processing.Indeed, a transformation onto a common quantity metric does not seem necessary to allow interactions between neural processing or perception of different quantities.A more nuanced interpretation of A Theory Of Magnitude is simply that there are related representations of different quantities, which many aspects of our results support.We thus suggest that behavioral interactions result from widespread interactions between distinct neural populations responding to different quantities within high-dimensional and distributed population codes for different quantities.

E
.Hendrikx et al.

Fig. 2 .
Fig. 2. Preferred quantity of voxels throughout the quantity map hierarchies.(A) Numerosity preferences for voxels with over 0.3 variance explained by the numerosity response model.Numerosity map borders are shown as black lines.Timing map borders are shown as dashed white lines.(B) Duration preferences for voxels with over 0.2 variance explained by the timing response model.Timing map borders are shown as white lines.Numerosity map borders are shown as dashed black lines.(C) As (B) but with period preferences.The light shaded region is outside the fMRI recording volume.Represented are two example participants, for all participants see Supplementary Figs.2-3.

Fig. 4 .
Fig.4.Consistent positive correlations coefficients between preferred numerosity and timing around IPS. (A) Correlation coefficients between preferred numerosity and preferred duration for each numerosity map.(B) As in (A), but for correlations between preferred numerosity and preferred period.The distributions of correlation coefficients were significantly above zero only in some numerosity maps around the IPS (NPO, NPCM, starred).(C) Correlation coefficients between preferred numerosity and preferred duration for each timing map.(D) As in (C), but for correlations between preferred numerosity and preferred period.The distributions of correlation coefficients were significantly above zero only in one timing map around the IPS (TPO).Brackets and stars show significant differences in multiple comparisons with maps grouped across hemisphere: all brackets to the left of the star were significantly different from all brackets to the right of the star (p < 0.05).Squares show the median correlation coefficient across all examples of the same map, error bars show the 95 % confidence interval of the median.Colored markers follow the formats of Fig.3.(For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Fig. 5 .
Fig. 5. Schematic overview of the angle computation in a numerosity map example.(A) For each voxel in the numerosity (black border) or timing (white border) maps, we calculated vectors of each quantity preference's progression (arrows).In overlapping areas, we calculated the angle between these vectors (gray shaded, circular area).(B) We averaged these angles over all overlapping voxels within a map.

Table 1
Descriptive and test statistics for testing the amount of overlap with any other map of the other quantity against a chance level overlap.interval of the median computed from 1000 bootstrap iterations.These are the outcomes of a one-sided Wilcoxon signed-rank test (FDR corrected for multiple comparisons).

Table 2
Descriptive and test statistics for testing the distribution of Spearman correlation coefficients between quantity preferences in all overlapping maps against 0.

Table 2 (
continued ) 95 % confidence interval of the median computed from 1000 bootstrap iterations.These are the outcomes of a two-sided Wilcoxon signed-rank test (FDR corrected for multiple comparisons).No statistical comparisons were done for TLS, as n ≤ 1.