Gradual change of cortical representations with growing visual expertise for synthetic shapes

An essential aspect of visual object recognition is the processing of visual shapes. The neural substrate of shape processing includes the ventral visual pathway, which in humans extends over the ventral occipitotemporal cortex from the occipital pole to the lateral occipital cortex, fusiform gyrus, and beyond (reviewed by Bi et al., 2016; Grill-Spector & Weiner, 2014; Kravitz et al., 2013; Weiner & Zilles, 2016). Functional imaging studies of ventral occipitotemporal cortex reveal intriguing functional anatomy, with responsiveness to specific object categories (e.g., faces, scenes, body parts) changing systematically over the cortical surface along several large-scale anatomical gradients (e.g., animate-inanimate, large-small, feature-whole, or perception-action; Freud et al., 2017; Grill-Spector & Weiner, 2014; Grill-Spector et al., 2004; Konkle & Oliva, 2012; Wurm & Caramazza, 2022; Yildirim et al., 2019).

Experience and learning improve object recognition performance, and also modify shape processing in the ventral occipitotemporal cortex. Indeed, functional imaging evidence shows that particular visual expertise—being able to identify and categorize visually similar objects of a particular kind—often entails moderate but anatomically distributed changes in the pre-existing responsiveness to shape (reviewed by Bukach et al., 2006; de Beeck & Baker, 2010; Gauthier & Tarr, 2016; Harel et al., 2013). This has been established by comparing novices and experts for identifying particular categories of natural objects (e.g., birds, mushrooms, minerals, degraded images; Cetron et al., 2019; Connolly et al., 2012; Duyck et al., 2021; Freud et al., 2017; Martens et al., 2018; McGugin et al., 2012; Roth & Zohary, 2015), as well as by comparing observers before and after they have learned to categorize initially unfamiliar synthetic shapes (e.g., computer-generated “greebles,” “spikies,” or “ziggerins”; Brants et al., 2011; de Beeck et al., 2006; Gauthier et al., 1999; A. C.-N. Wong et al., 2009; Y. K. Wong et al., 2012; Yue et al., 2006).

Here, we map the cortical representation of synthetic visual objects and track gradual changes as initially unfamiliar objects become progressively familiar with learning. We wondered how pre-existing shape representations would accommodate and integrate novel synthetic objects. We further wondered whether representational changes would be specific to learned objects or extend also to other objects of the same kind. To explore these questions, we analyzed “representational similarity” of spatiotemporal BOLD patterns (Haxby, 2012; Kriegeskorte, Mur, Ruff, et al., 2008), which offers a potentially sensitive measure for the information encoded in neural activity and may also be related to similarity as perceived by human observers (Charest & Kriegeskorte, 2015; Collins & Behrmann, 2020; Nestor et al., 2016).

Most previous studies of visual expertise identified cortical sites associated with a particular object category by comparing BOLD activity either between novices and experts or before and after learning. We extend this work in three ways: firstly, by establishing representational distance at the level of object exemplars rather than object categories; secondly, by monitoring gradual changes as observers gain familiarity with object exemplars; and thirdly, by analyzing changes in the diversity of multivariate BOLD activity. Few previous studies have attempted to resolve shape representations in such detail (Brants et al., 2016; Duyck et al., 2021; Eger et al., 2008; Visconti di Oleggio Castello et al., 2021). To progress fine-grained analysis of representational geometry, we developed synthetic shapes for which visual expertise is acquired comparatively slowly (Kakaei et al., 2021) and took advantage of a numerically tractable method for linear discriminant analysis in O(⁠$103$⁠)-dimensional multivariate activity (DLDA; Yu & Yang, 2001).

Our results showed view-invariant representations of shape over surprisingly extensive regions of the ventral occipitotemporal cortex, including the fusiform gyrus, lateral occipital areas, and primary visual cortex. Representational distances were high from the start, even before learning, suggesting that new visual expertise was accommodated and encoded within pre-existing representations. However, shapes that appeared repeatedly (and were memorized by observers) and shapes that appeared just once (and were ignored) diverged dramatically, in terms of their cortical representations, while visual expertise was being acquired and consolidated.

Eight healthy observers (4 female and 4 male; aged 25 to 32 years) took part in behavioral training (“sham experiment,” one session per observer), the functional imaging experiment (“main experiment,” six scanning sessions per observer), and a final behavioral assessment (two sessions). All observers were paid and gave informed consent. Ethical approval was granted under Chiffre 30/21 by the ethics committee of the Faculty of Medicine of the Otto-von-Guericke University, Magdeburg.

In both sham and main experiments, observers viewed sequences of 200 recurring and non-recurring objects (see below and Fig. 1A) and attempted to classify each object as “familiar” or “novel” (by pressing the appropriate button). Over the course of multiple sessions, observers gradually became familiar with recurring objects and thus became able to distinguish them from non-recurring objects. Objects of the sham experiment were two-dimensional shapes, whereas objects of the main experiment were rotating, three-dimensional shapes (see below and Fig. 1A).

The main experiment extended over 3 successive weeks, with three sessions on separate days of both the 1st and 3rd week (no sessions took place in the 2nd week). The experiments of the 1st and 3rd week differed in four aspects: sequence type (structured or unstructured), the set of recurring objects, object color (red or blue), and responding hand (left or right). All aspects were counterbalanced across observers.

After the three scanning sessions of a week, observers participated in an additional behavioral session to confirm that they had in fact become familiar with every recurring object. Specifically, they performed a spatial search task in which they pointed out recurring target objects among non-recurring distractor objects (Kakaei et al., 2021). In addition, observers were offered the opportunity to voice anything they might have noticed about the experiment.

Complex three-dimensional objects were computer-generated and presented as described previously (Kakaei et al., 2021). A movie can be viewed under this LINK. All objects were highly characteristic and dissimilar from each other as confirmed computationally in terms of vector distances between depth maps (Kakaei et al., 2021). Objects were presented every $3s$⁠, with $2.5s$ viewing and $0.5s$ transition time (Fig. 1A). Objects were shown from all sides and, after appearing at an arbitrary angle, revolved smoothly for one full turn (period $2.5s$⁠, frequency $0.4Hz$⁠, angular frequency $144∘​/​s$⁠) about one of several axes in the frontal plane (⁠$−45∘$⁠, $0∘$⁠, $45∘$⁠, clockwise or counter-clockwise). Axes and directions were counterbalanced for each object, and initial viewing angles were chosen randomly (Fig. 1B). All stimuli were generated with MATLAB (The MathWorks, Inc.), presented with the psychophysics toolbox (Brainard, 1997), and viewed in a mirror mounted to the MR head coil (screen resolution $960×720$ pixels, frame rate $60Hz$⁠, subtending approximately $8∘​×6∘$ of visual angle, average luminance $50Cd​/​m2$⁠, background luminance $5Cd​/​m2$⁠). Observers responded with the right or left index finger on an MR-safe response box.

Fifteen objects recurred many times during three sessions (“recurring” objects), whereas other objects appeared exactly once (“non-recurring” or “singular” objects). As mentioned, observers classified every object as either “familiar” or “unfamiliar” by pressing a button during its presentation. Over the course of three sessions, all observers gradually became familiar with the “recurring objects” (see below). The average time-course of learning, as established by a simplified signal detection and reaction-time (RT) analysis, is shown in Figure 1C.

Every session comprised six sequences (“runs”), each lasting $600s$ and presenting $180$ “recurring” and $20$ “non-recurring” objects (⁠$200$ objects in total). As there were 15 different recurrent objects, each such object was seen $12±1.9$ times during every sequence. Over the three sessions (or 18 sequences), each recurring object appeared at least $190$ times each (mean $±$ S.D.: $216±9$⁠), whereas non-recurring objects appeared only once. Altogether, there were $3,240$ presentations of recurring objects (⁠$3×6×180$⁠) and $360$ presentations of non-recurring objects (⁠$3×6×180$⁠).

Presentation sequences started with a random recurring object and continued randomly to one of the possible next objects, with neither immediate repetitions (⁠$X→X$⁠) nor direct returns (⁠$X→Y→X$⁠) being allowed. Sequences comprised $200$ objects, of which $180$ were recurring and $20$ objects non-recurring and were interspersed at random intervals. Object sequences were post-selected such as to counterbalance the number of appearances of every recurring object in every session.

All observers performed the experiment twice in the scanner, once during the 1st week and again during the 3rd week of the main experiment (so that 8 observers provided 16 data sets). As mentioned, the 2 weeks differed in terms of the recurring objects and the presentation sequence. “Structured” sequences exhibited predictive sequential dependencies (3 possible recurring next objects), whereas “unstructured” sequences did not (14 possible recurring next objects, see Kakaei et al., 2021 for details). As a result, the repetition latency (i.e., the latency of successive presentations of the same object) was $5.5±15$ (median and S.D.) for “structured” and $10.5±11$ for “unstructured” sequences. Further aspects and effects of sequence structure are reported and discussed in detail in a companion paper.

To verify that recurring objects had become familiar to observers, every observer performed 60 trials of a spatial search task with 3 recurring and 9 non-recurring objects. The 12 objects were positioned randomly in a $3×4$ array and were presented for 30 s while rotating in three dimensions (as in the main experiment). After each presentation, observers indicated the recurring object positions with the computer mouse. Performance was consistently above 95% correct.

All magnetic-resonance images were acquired on a 3T Siemens Prisma scanner with a 64-channel head coil. Structural images were T1-weighted sequences (MPRAGE TR = 2,500 ms, TE = 2.82 ms, TI = 1,100 ms, $7∘$ flip angle, isotropic resolution $1×1×1mm$ and matrix size of $256×256×192$⁠). Functional images were T2*-weighted sequences (TR = 1,000 ms, TE = 30 ms, $65∘$ flip angle, resolution of $3×3×3.6mm$ and matrix size of $72×72×36$⁠). Field maps were obtained by gradient dual-echo sequences (TR = 720 ms, TE1 = 4.92 ms, TE2 = 7.38 ms, resolution of $1.594×1.594×2mm$ and matrix size of $138×138×72$⁠).

Our approach to fMRI analysis was influenced by recent advances in comparing uni- and multivariate responses of corresponding voxels between different observers (e.g., Kumar et al., 2022; Nastase et al., 2019). The local correlation structure of voxel response, which is similar in different observers, provided the basis for our functional parcellation (Dornas & Braun, 2018). The parcellation obviated “searchlight” strategies by defining for all observers corresponding brain “parcels” with corresponding episodes of high-dimensional (O(1000)) multivariate activity.

The fMRI pre-processing procedure was similar to that published previously (Dornas & Braun, 2018). First, DICOM files were converted into NIFTI format using MRICRON (MRICRON Toolbox, Maryland, USA, NIH). Then, brain tissues were extracted and segmented using BET (Smith, 2002) and FAST (Zhang et al., 2001). Field map correction, head motion correction, spatial smoothing, high-pass temporal filtering, and registration to structural and standard images were performed with the MELODIC package of FSL (Beckmann & Smith, 2004).

Field map correction and registration to structural image were carried out using Boundary-Based Registration (BBR; Greve & Fischl, 2009). MELODIC uses MCFLIRT (Jenkinson et al., 2002) to correct for head motion. Spatial smoothing was performed with SUSAN (Smith & Brady, 1997), with full width at half maximum set at FWHM $=5mm$⁠. To remove low-frequency artifacts, we applied a high-pass filter of the cut-off frequency $f=0.01$Hz, that is, oscillations/events with periods of more than $100$ s were removed. To register the structural image to Montreal MNI152 standard space with isotropic $2mm$ voxel size, we used FLIRT (FMRIB’s Linear Image Registration Tool; Jenkinson & Smith, 2001; Jenkinson et al., 2002) with 12 degrees of freedom (DOF) and FNIRT (FMRIB’s Nonlinear Image Registration Tool) to apply the non-linear registration. To further reduce artifacts arising from head motion, we applied despiking with a threshold of $λ=100$ using BrainWavelet toolbox (Patel et al., 2014). Later, we regressed out the mean CSF activity as well as 12 DOF translation and rotation factors predicted by a motion correction algorithm (MCFLIRT). Afterward, the time series of each voxel was detrended linearly and whitened (with Matlab functions “detrend” and “zscore”).

Finally, the $160,099$ voxels of MNI152 space were grouped into $758$ functional parcels according to the MD758 atlas (Dornas & Braun, 2018). Each functional parcel is associated with an anatomically labeled region of the AAL atlas (Tzourio-Mazoyer et al., 2002) and comprises approximately $200$ voxels or approximately $1.7cm3$ of gray matter volume (⁠$212±70$ voxels, range $45$ to $462$ voxels). Parcels were defined for a small population of observers such as to maximize signal covariance within and minimize covariance between parcels in the resting state. In contrast to other parcellation schemes, this was based exclusively on the (typically strong) functional correlations within each anatomical region and disregarded the (typically weak) correlations between different anatomical regions. The MD758 parcellation offers superior cluster quality, correlational structure, sparseness, and consistency with fiber tracking, compared to other parcellation schemes of similar resolution (Albers et al., 2021; Dornas & Braun, 2018).

To study the neural representation of objects, we extracted the multivoxel activity pattern at $Nt=9$ time points following object onset. In a functional parcel with $Nvox$ voxels, this response pattern constituted a point (or vector) in an $Ndim$-dimensional space, where $Ndim=Nt⋅Nvox$ (Fig. 2A). To identify parcels with significant selectivity for individual recurring objects, we employed a representational similarity analysis (RSA; Kriegeskorte, Mur, & Bandettini, 2008) (Fig. 2B). This analysis uses the standardized Euclidean (Mahalanobis) distance between responses in a high-dimensional space to examine the separability of neural object representations as a function of learning, or object type (recurring or non-recurring), or both. Over all 758 parcels, response dimensionality was $Ndim=1,911±634$ (mean and standard-deviation), with a range from $405$ (Calcarine-L 329, with $45$ voxels) to $4,113$ (Postcentral-R-484, with $457$ voxels).

Our approach to RSA differed from previous work in some respects. Firstly, we analyzed high-dimensional spatiotemporal patterns of BOLD activity (⁠$200$ voxels $×$$9$ s, or O(⁠$103$⁠) dimensions) in non-overlapping gray matter volumes (⁠$758$ functional subdivisions of $90$ anatomical regions, averaging $1.7 cm3$⁠; Dornas & Braun, 2018). Other studies have used lower-dimensional spatial activity patterns in overlapping searchlight volumes (⁠$O(102)$ voxels or dimensions, covering $0.25$ to $1.0cm3$⁠; Kriegeskorte et al., 2006). Secondly, we employed multi-class linear discriminant analysis (“direct linear discriminant analysis,” DLDA; Yu & Yang, 2001), rather than pairwise discriminability or one-versus-all discriminability (e.g., Hung et al., 2005; Liu et al., 2009). With these modifications, RSA revealed representational geometry at the level of object exemplars, as well as gradual changes in this geometry over sessions and runs.

To analyze the response variance that discriminates $κ=15$ recurring objects, at most $(κ−1)$-dimensions are required. Restricting the analysis to $14$ principal components of the response could potentially have neglected smaller but more discriminating components. Accordingly, we performed a Linear Discriminant Analysis (LDA), which amounts to a “supervised” principal component analysis (PCA) and yields the $(κ−1)$-dimensional orthonormal subspace $S$ that optimally discriminates the $κ$ response classes. Here, optimality is defined as simultaneously minimizing within-class variance and maximizing between-class variance of responses.

The results of LDA and PCA showed considerable commonality. Over the 758 parcels, the first 14 principal components captured $53±7%$ (mean and S.D.) of the total response variance, whereas the 14-dimensional subspaces $S$ captured $33±7%$ of the total variance (or $61±6%$ of the principal component variance). Almost all of the subspace variance overlapped with the principal component variance (i.e., $88±5%$ of subspace variance projected into the space of the first 14 principal components, while the remaining $12±5%$ projected into the space of the remaining principal components).

Similar numbers were obtained for the 124 identity-selective parcels. The first 14 principal components captured $57±6%$ (mean and S.D.) of the total response variance, and subspaces $S$ captured $38±6%$ of the total variance (or $67±4%$ of the principal component variance). Almost all of the subspace variance (⁠$91±3%$⁠) overlapped with the first 14 principal components. In summary, Linear Discriminant Analysis captured the useful (discriminating) part of correlated variance and distributed this variance more uniformly over its 14 orthonormal dimensions (⁠$6±3%$ per dimension) than principal component analysis could (⁠$4±6%$ per dimension).

A numerically tractable procedure for identifying the optimal subspace $S$ is available in terms of “direct LDA” or DLDA (Ye et al., 2006; Yu & Yang, 2001). Briefly, this method first diagonalizes between-class variance to identify $κ−1$ discriminative eigenvectors with non-zero eigenvalues, next diagonalizes within-class variance, and finally yields a rectangular matrix for projecting activity patterns from the original activity space (dimensionality $Ndim$⁠) to the maximally discriminative subspace $S$ and back. As this method is linear and relies on all available degrees of freedom, its results are deterministic. An important feature of this particular algorithm is that within-class variance is maintained near unity for all classes, by means of a suitable scaling of the subspace dimensions. The link github.com/cognitive-biology/DLDA provides a Matlab implementation of DLDA.

Activity patterns $xjk$ associated with trials $k$ were analyzed in the maximally discriminative subspace $S$⁠. The normalized amplitude $ak=1κ−1∑j=1κ−1xjk2$ of such patterns exhibited an average value of $〈a〉=0.99$⁠. The normalized distance $dkl=1κ−1∑j=1κ−1(xjk−xjl)2$ between patterns associated with trials $k$ and $l$ measured on average $〈d〉=1.40$⁠, consistent with distance expected between random patterns of this amplitude (⁠$2$⁠). Averaging over trials $k$ produced normalized response amplitudes $A=〈ak〉k$⁠. Averaging over pairs of trials $k$⁠, $l$ separated by a given latency $l−k$⁠, produced normalized response distances $D(l−k)=〈dkl〉k,l$⁠.

The patterns from successive trials exhibited a weak temporal correlation, with approximately $5%$ smaller distances at delays below $4$ trials and approximately $2%$ larger distances at delays ranging from $6$ to $15$ trials (see Supplementary Fig. S1A, B). Comparing pairs of trials with different types of objects, we observed approximately $3%$larger response distances $D$ (at all latencies) for the same recurring objects than for either different recurring or non-recurring objects (Supplementary Fig. S1C). Differential response amplitudes $A$ increased marginally with latency, because response amplitudes tended to increase slightly over the course of each run (Supplementary Fig. S1D). This trend was evident for all types of objects and with both “structured” and “unstructured” sequences. In other words, the effect of object type on multivariate hemodynamic responses was limited to response distances and did not extend to response amplitudes. Thus, our data provided no evidence for “repetition suppression.”

For certain analyses (Sections 2.5.8 and 2.5.9), we established for each parcel $w$ the average delay-dependent distance $Tw(Δk)=〈dw,u,r(Δk)〉u,r$ between patterns with a relative delay of $Δk$ trials, where the average was taken over subjects $u$ and runs $r$⁠. The time-course $Tw$ allowed us to discount temporal correlations by computing $dw,u,rcorrected(Δk)=dw,u,r(Δk)−Tw(Δk)+〈Tw(Δk)〉Δk$⁠, where $〈Tw(Δk)〉Δk$ is the average value over delays $Δk$⁠.

Our observations comprised approximately $200$ activity patterns for each of the $15$ recurring object classes (per observer and condition). To allow for cross-validation, we randomly divided these patterns in a larger “training set” (⁠$90%$ or $190±7.7$ per object class) and a smaller test set (⁠$10%$ or $22±0.9$ per object class) (Fig. 2B). Note that the “training set” comprised exclusively activity patterns associated with recurring objects. To reduce the variability introduced by random test sets, this selection was repeated $Nr=20$ times and all statistical measures described below represent the average over repetitions. As illustrated in Figure 2C, in the discriminative subspace $S$⁠, we compared the $ni$ test set exemplars $xki$ (where $k=1,…,ni$⁠) of class $i$ to the centroids $cjtrain$ established for the training exemplars of class $j$⁠. To compute Mahalanobis distances and variance ratios (see below), we compared test set exemplars $xki$ of class $i$ to the centroids $cjtest$ of test set exemplars of class $j$⁠.

We used three measures for this comparison, all with comparable results. Firstly, the nearest class centroid $citrain$ to each pattern exemplar $xki$ was identified to establish a matrix of classification probabilities $P(j​|i)$ (probability that an exemplar of class $i$ is nearest to the centroid of class $j$⁠), also known as “confusion matrix,” as well as the “classification accuracy” $α=∑iP(i​|i)P(i)$⁠, which is the probability that the nearest centroid is the correct one.

Secondly, for each pair of object classes $(i,j)$⁠, object exemplars $xki$ and $xkj$ from the test set were projected onto the line connecting the two test set centroids, $citest$ and $cjtest$⁠, and a pairwise discriminability/dissimilarity/Mahalanobis distance $δi,j$ was computed from the means, $μi$ and $μj$⁠, and variances, $σi2$ and $σj2$⁠, of these projections, as $δi,j=|μi−μj|12(σi2+σj2)$⁠. The average over all pairs of object classes was computed as $δ=2κ(κ−1)∑i,jδi,j$⁠.

Thirdly, given class centroids $citest$ and overall centroid $ctest$⁠, we computed the Euclidean distances $dki=‖xki−citest‖$ between exemplars $xki$ and class centroid $citest$ and, for each object class $i$⁠, the “sum of squares” as $SSWi=∑k=1nidki2$⁠. The “within-class” variance of all classes was computed as $SSW=1N∑i=1κSSWi$⁠, where $N=∑i=1κni$⁠. Similarly, from the Euclidean distances $di​=‖citest−ctest‖$ between individual and overall centroids, we computed “between-class” variance $SSB​=1N∑i=1κnidi2$⁠. From the Euclidean distances $dki​=‖xkitest−ctest‖$ between exemplars and overall centroid, we computed “total” variance $SST​=1N∑i=1κ∑k=1nidki2$⁠. Variances $SSW$⁠, $SSB$⁠, and $SST$ are also denoted, respectively, $SSsame$⁠, $SSdiff$⁠, and $SSfam$ further below. To quantify the discriminability of classes, the variance ratio $Fidentity​=SSB(N−κ)​/​SSW(κ−1)$ provided a non-parametric multivariate statistic (PERMANOVA; Anderson, 2001). The average within-class and between-class dispersion per dimension could be estimated as $σW​=SSW​/​(N−κ)$ and $σB=SSB​​/​(κ−1)$⁠, respectively.

To test for statistical significance, we computed average classification performance (in terms of both classification accuracy $αobs$ and f-ratio $Fobs$⁠) over $Nr$ test sets, as well as over $103$ first-level permutations of object identities (in each of the $Nr$ test sets). In principle, we could have tested an “individual null” hypothesis for every parcel and every data set, namely, the probability of obtaining the observed performance $αobs$ (or $Fobs$⁠) purely by chance. Instead, we computed the “minimum statistic” $m=minkαk$ (or $m=minkFk$⁠) over data sets $k$⁠, as well as over $105$ second-level permutations (drawn randomly from the first level permutations) and tested the “global null” hypothesis, namely, the probability $pn(m)$ of obtaining the observed minimum performance over $n$ data sets purely by chance. This computation was performed separately for each of the 2 conditions (8 data sets from 8 observers per condition) as well as for the union of conditions (16 data sets from 8 observers). When the “global null” hypothesis could be rejected, we inferred statistically significant classification performance in at least some data sets. Our threshold for significance was $pn⋆(m)<0.05$ after correction for multiple comparisons (⁠$758$ parcels and $2$ conditions) (Allefeld et al., 2016).

To summarize the results from all observers and conditions, we used a “prevalence analysis” (Allefeld et al., 2016). Prevalence $γtrue$ is the fraction of significant performance over $n=16$ data sets. To test the “prevalence null” hypothesis that $γtrue$ is below a threshold $γ0​=0.5$⁠, an upper bound for $P(γtrue​<γ0)$ was obtained from the probability $pn⋆(m)$ of the minimum statistic over $n=16$ data sets, after correction for multiple comparisons:

This was the criterion used to label parcels as “identity selective.” Threshold prevalence $γ≃0.5$ corresponded to corrected probability $pn⋆(m)≃0.0012$ and minimal accuracy of $6.67%$ (i.e., near chance).

Additionally, we computed $γest$ as the largest value for which the “prevalence null” hypothesis could be rejected from

where $pn⋆(m)$ is the corrected minimum probability, $n=16$ the number of data sets, and $α=0.05$ the significance threshold.

Although recurring and non-recurring objects were comparable and generated in the same way, it seemed possible that neural representations might discriminate the class of 15 recurring objects from the class of 360 non-recurring objects. Indeed, the two classes became discriminable after observers had learned to classify recurring objects as “familiar” and non-recurring objects as “novel.” Accordingly, we considered this discriminability a representation of “novelty.”

To assess the neural representation of “novelty,” we divided non-recurring and recurring objects into two sets of unequal size (approximately $N=216×15$ recurrent or “familiar” exemplars vs. $M=360$ non-recurrent or “novel” exemplars). From the Euclidean distances $dk​=‖xk−c‖$ between test set exemplars $xk$ and centroids $cfam​=1N∑k=1Nxk$ or $cnov​=1M∑k=1Mxk$⁠, we obtained “within-class” variance $SSW​=SSfam​+SSnov$⁠, where $SSfam​=1N+M∑k=1Ndk,fam2$ and $SSnov​=1N+M∑k=1Mdk,nov2$⁠. From distances $dfam​=‖cfam​−ctot‖$ and $dnov​=‖cnov​−ctot‖$ between class centroids and overall centroid $ctot​=NN+Mcfam​+MN+Mcnov$⁠, we obtained “between-class” variance $SSB​=SSnovfam=NN+Mdfam2+MN+Mdnov2=$$NM(N+M)2(cfam​−cnov)2$⁠. Finally, from distances $dk=$$​‖xk​−ctot‖$ between exemplars and overall centroid, we obtained total variance $SST=1N+M∑k=1N+Mdk2$⁠. To quantify the discriminability of non-recurring and recurring objects, we formed the variance ratio $Fnovelty​=SSB(N+M​−2)​/​SSW$ (Anderson, 2001). Average within-class and between-class dispersion per dimension was obtained from $σW​=SSW​/​(N+M−2)$ and $σB=SSB$⁠, respectively.

To assess changes in neural representations over the course of the experiment, while also allowing for cross-validation, we divided all recurring object presentations into five successive “batches” $B1,B2,…$⁠, each with $20%$ of the presentations (Fig. 2D). In this way, we could select “test sets” for cross-validated DLDA from one particular batch, while retaining all other presentations as a “training set.” As every recurrent object was presented $210±9$ times over all sessions, a batch would comprise $42±1.8$ presentations, a test set $21±0.9$⁠, and a training set $189±8.1$ presentations. To reduce the variance deriving from test set selection, we repeated the random selection $Nr​=20$ times and averaged over repetitions.

To quantify representational changes over the course of learning, we computed the variance ratios $Fm,w,uidentity$ for each temporal window or batch $m$⁠, identity-selective parcel $w$⁠, and data sets $u∈{1,…,16}$⁠. We formed the average ratio over 16 data sets, $Fm,widentity​=〈Fm,w,uidentity〉u$⁠, and assessed statistical significance by shuffling (⁠$103$ permutations) the identity of recurring objects to obtain the distribution of variance ratios due to chance or data structure. The mean $μm,w$ and variance $σm,w2$ of this distribution could also be used to convert $Fm,widentity$ into z-score values $Zm,widentity=(Fm,widentity​−μm,w)/σm,w$⁠.

Additionally, we performed a regression analysis and quantified representational changes in terms of linear trends. Specifically, we determined a “rate” parameter $βwidentity$ by fitting a linear mixed-model $Fm,w,uidentity=$$​β0,w​+βwidentitym+ξ0,w,u​+ξ1,w,um+ϵm,w,u$ with data sets $u$ as the grouping variable, where $β0,w$ was a fixed-effect coefficient, $ξ0,w,u$ and $ξ1,w,u$ were random effect coefficients, and $ϵm,w,u$ was residual error.

Similarly, to assess whether neural representations of non-recurring objects change with learning, we divided all object presentations (recurring and non-recurring) into five successive “batches” $B1,B2,...$⁠, each with $20%$ of the presentations (Fig. 2D), to obtain variance ratios $Fm,w,unovelty$ for each temporal window or batch $m$⁠, identity-selective parcel $w$⁠, and data sets $u∈{1,…,16}$⁠. After averaging over 16 data sets, $Fm,wnovelty=〈Fm,w,unovelty〉u$⁠, we assessed statistical significance by shuffling (⁠$103$ permutations) the identity of recurring and non-recurring objects to obtain the distribution of variance ratios due to chance or data structure. The mean $μm,w$ and variance $σm,w2$ of this distribution were used to convert $Fm,wnovelty$ into z-score values $Zm,wnovelty​=(Fm,wnovelty​−μm,w)/σm,w$⁠.

Additionally, we performed a regression analysis to establish linear trends. Changes in the representation of object “novelty” were assessed by fitting the “rate” parameter $βwnovelty$ in a linear mixed-model $Fm,w,unovelty​=β0,w​+βwnoveltym+ξ0,w,u​+ξ1,w,um+ϵm,w,u$⁠, with data sets $u$ as the grouping variable, where $β0,w$ was a fixed-effect coefficient, $ξ0,w,u$ and $ξ1,w,u$ were random effect coefficients, and $ϵm,w,u$ was a residual error.

To establish linear trends $Fm=〈Fm,w,u〉w,u$ (of either identity and novelty) that average over both parcels $w$ and data sets $u$⁠, we obtained a rate parameter $β1$ by fitting linear mixed-model $Fm,w,u​=β0​+β1m+ξ0,w,u​+ξ1,w,um+ϵm,w,u$ with both parcels and data sets as grouping variables.

In the cross-validated analyses described above, subspaces $S$ differed slightly between different batches (and training sets). To analyze the geometry of neural representations in a stable framework, we repeated some analyses in fixed subspaces $S$ that reflected all observations (i.e., all recurring activity patterns $xk$⁠). In the fixed subspace, we calculated the normalized amplitude $ak​=‖xk‖/​κ−1=∑j=1κ−1xjk2/​κ−1$ of individual patterns $k$ and the normalized pairwise distance $dkl​=‖xk​−xl‖/​κ−1=∑j=1κ−1(xjk​−xjl)2/​κ−1$ between two patterns $k$ and $l$⁠.

For each parcel $w$⁠, data set $u$⁠, and run $r$⁠, we obtained the average amplitude $Aw,u,rtot​=1N+M∑k=1N+Mak$ of all patterns, the average amplitude $Aw,u,rfam​=1N∑k=1Nak$ of recurring patterns, and the average amplitude $Aw,u,rnov​=1M∑k=1Mak$ of non-recurring patterns. Similarly, we obtained the average pairwise distance $Dw,u,rtot​=2(N+M)(N+M−1)$$∑k=1N+M∑l=kN+Mdkl$ between all patterns, the average distance $Dw,u,rnov=2M(M−1)∑k=1M∑l=kMdkl$ between non-recurring patterns, the average distance $Dw,u,rfam​=$$2N(N−1)∑k=1N∑l=kNdkl$ between recurring patterns, and the average distance $Dw,u,rnovfam​=1MN∑k=1M∑l=1Ndkl$ between pairs comprising one recurring and one non-recurring pattern. For recurring patterns, we further obtained the average distance $Dw,u,rsame=2N(N​/​κ−1)∑i=1κ∑k=1ni∑l=knidkl$ between pairs of recurring patterns in the same class and the average distance $Dw,u,rdiff​=1N(N−N​/​κ)∑i=1κ∑k=1ni$$∑l=1N−nidkl$ between pairs in different classes. All distances were corrected for the temporal auto-correlation by subtracting the time course of $Tw(i,j)$⁠, as described above.

As described further above, the distances between individual activity patterns and different centroids—such as $ctot$⁠, $cnov$⁠, and $cfam$—yielded total variance $SST​=SStot$⁠, within-class variance $SSW​=SSfam​+SSnov$⁠, and between-class variance $SSB​=SSnovfam$⁠. For recurring patterns, distances to individual class centroids $ci$ and overall centroid $cfam$ yielded total variance $SST​=SSfam$⁠, within-class variance $SSW​=SSsame$⁠, and between-class variance $SSB​=SSdiff$⁠.

These values were computed for each parcel $w$⁠, observer $u$⁠, and run $r$⁠, in order to obtain variance fractions $Fw,u,rfam​=SSfam​/​SStot$⁠, $Fw,u,rnov​=SSnov​/​SStot$⁠, $Fw,u,rnovfam=$$​SSnovfam​/​SStot$⁠, $Fw,u,rsame​=SSsame​/​SSfam$⁠, and $Fw,u,rdiff=SSdiff/$$​​SSfam$⁠, as well as variance ratios $Rw,u,ridentity​​=SSdiff(N−κ)​/$$SSsame(κ−1)$ and $Rw,u,rnovelty​=SSnovfam(N+M−2)​/​(SSnov​​=$$SSfam)$⁠.

Fixed subspaces permitted us to assess representational changes between successive “runs.” To this end, we computed average amplitudes $Aw,u,r$⁠, distances $Dw,u,r$⁠, variances $SSw,u,r$⁠, and variance ratios $Fw,u,r$⁠, as described above, for each parcel $w$⁠, data set $u∈{1,…,16}$⁠, and run $r$⁠. Within each session $s$⁠, we assessed the changes of these parameters $Y​∈{A,D,SS,F}$ over runs $r′∈s$ by determining a “rate” parameter $βs$ for identity-selective $w$ and non-selective parcels $w′$⁠. Each $βs$ coefficient was acquired from a linear mixed-model $Yr′,w,u​=β0,s​​+βsr′+ξ0,w,u​+ξ1,w,ur′​+ϵr′,w,u$ with observers and parcels as grouping variables, where $β0,s$ was a fixed-effect coefficient, $ξ0,w,u$ and $ξ1,w,u$ were random effect coefficients, and $ϵ$ was residual error. The same approach was used to assess gradual changes over runs in the centroid-to-centroid distances $Dsame(r)$⁠, $ΔDsame(r)$⁠, $Dnov(r)$⁠, and $ΔDnov(r)$⁠. This served to test the statistical significance of linear rates $βs$ in each session. Sessions with significant rates are marked by stars in Figure 6.

We also assessed the stability of the representation of the 16 response classes (15 recurring and 1 non-recurring) over the course of the experiment. To this end, we compared the average representation in individual runs $r$ (centroids $Cr$ of responses to exemplars) to the average representation over all runs (centroids $Cave$⁠). For observer $u$⁠, identity-selective parcel $w$⁠, and object class $i$⁠, we calculated the Euclidean distance $Du,w,i,r$ between the relevant $Cr$ and $Cave$⁠, and also the difference $ΔDu,w,i,r$ between the relevant centroids from successive runs, $Cr$ and $Cr+1$⁠. After averaging over observers $u$⁠, identity-selective parcels $w$⁠, and object classes $i$⁠, we obtained $Dsame(r)$ and $ΔDsame$ for recurring objects and by $Dnov(r)$ and $ΔDnov(r)$ for non-recurring objects.

As a baseline for comparison, we also computed the distances $Du,w,i,r$ and differences $ΔDu,w,i,r$ that may be expected purely on the basis of response variance. To this end, we permuted the sequence of all $3,600$ trials, separately within each of the 16 response classes (15 recurring and 1 non-recurring) such as to obtain 18 “pseudo-runs” with $200$ trials each. Expectation values were obtained by repeating this $Nr​=1,000$ times.

We note that, in an $n$-dimensional hypersphere of unit radius, the average Euclidean distance between two random points is

with $dave≈1.4017$ for $n=14$⁠.

To visualize representational geometry in two dimensions, we randomly sampled $50$ response patterns to each of the recurring and non-recurring objects within the first and the last sessions and calculated a $1,600×1,600$ pair-wise distance matrix (⁠$Dw,u$⁠) for each identity-selective parcel $w$ and subject $u$⁠. We did not wish to average distance matrices over observers, as we did not expect the activity patterns of different observers to be comparable. To sidestep this difficulty, we permuted the order of recurring objects $100$ times and for each subject obtained an average matrix $D¯$ over permutations, which was then averaged over subjects. To visualize the representational geometry of identity in the first and the last session, we used multidimensional scaling (Matlab function mdscale, metric stress) to map the distances matrices for recurring objects (⁠$50$ exemplars from the first session and $50$ exemplars from the last session) into a two-dimensional space. To visualize the representational geometry of novelty, we restricted the distance matrix to non-recurring objects (⁠$50$ exemplars from the first session and $50$ exemplars from the last session) and just 3 of the 15 recurring objects (⁠$20$ exemplars from either session).

Observers viewed sequences of computer-generated objects, with each object shown for $2.5s$ while rotating in three dimensions (Fig. 1A, B, a movie may be viewed HERE). Over three sessions, observers viewed 3,600 objects in total, of which 3,240 were presentations of recurring objects (15 different objects, each appearing approximately 216 times) and 360 were presentations of non-recurring objects (360 different objects, each appearing once). The display was intended to be sufficiently intriguing to remain interesting over 3 successive days. To this end, presentations never repeated exactly. Observers were required to classify each object as “familiar” (recurring) or “novel” (non-recurring). The task performance improved as observers became increasingly familiar with recurring objects, as illustrated in Figure 1C. Over the first 600 presentations, classification performance improved approximately from $50%$ correct (chance) to $85%$ correct, and reaction times decreased approximately from $1.65s$ to $1.25s$⁠. Over the remaining 3,000 presentations, performance improved further to approximately $90%$ correct and reaction times decreased further to approximately $0.95s$⁠. After three sessions, all observers were “familiar” with all recurring objects and could pick them out from an array of distractor objects.

All sessions were performed in an MRI scanner while whole-brain functional imaging data were being collected. In the following, we report the results of three types of analyses. First, we describe the cortical areas in which multivariate BOLD activity encodes information about the identity of recurring objects (“object identity”), as determined by cross-validated analyses of entire data sets (3 sessions per observer). Second, we describe changes in cortical representations over coarse time intervals, by means of cross-validated analyses of successive parts of the data sets (3 sessions divided into 5 batches). These changes pertain to the encoding of both recurring objects and the distinction between recurring and non-recurring objects (“object novelty”). Third, we describe changes in representations over finer time intervals (3 sessions divided into 18 runs), by foregoing cross-validation and adopting a fixed reference frame. These finer intervals confirm the results from coarse intervals but reveal more details about the geometry of neural representations and their development over time.

To assess the extent to which multivariate neural responses to recurring objects encoded object identity, we relied on optimal linear classifiers combined with cross-validation (“direct linear discriminant analysis,” DLDA, see Methods for details). Specifically, we quantified the “identity” information in multivariate responses of every parcel $w∈{1,…,758}$ and data set $u∈{1,…,16}$ in terms of classification accuracy $αw,u$⁠, average pairwise dissimilarity $δw,u$⁠, and the ratio of between-class and within-class variance $Fw,u$⁠. All three measures proved highly correlated and supported similar conclusions. For example, Figure 3B illustrates the correlation of classification accuracy $αw,u$ and variance ratio $Fw,u$ (⁠$ρ=0.94$⁠, $p<0.001$⁠). The correlations of $aw,u$ and $δw,u$ (⁠$ρ=0.95$⁠, $p<0.001$⁠), and of $δw,u$ and $Fw,u$ (⁠$ρ=0.98$⁠, $p<0.001$⁠) were comparably strong. The results of individual observers from the two experimental conditions (structured and unstructured object sequences) were highly similar as well, demonstrating test-retest consistency (Supplementary Fig. S2).

For most parcels, the results from different observers showed considerable variability. Whereas a few parcels exhibited significant accuracy $αw,u$ and variance ratio $Fw,u$ in all data sets (e.g., Calcarine 331), in many parcels the representation of object identity was significant only in some data sets (e.g., Parahippocampus 325) (Fig. 3B). Global significance was assessed by comparing the minimal accuracy or variance ratio over the 8 data sets from one condition (structured or unstructured) to the minimal values obtained with shuffled data (red ellipse in Fig. 3C, see Methods for details).

Minimal classification accuracy $αw$ was significant in 17% of all parcels (128 of 758 parcels) in the structured sequence condition and in 19% of parcels (146 of 748) in the unstructured condition (⁠$p⋆0.05$⁠, corrected for multiple comparisons), when compared to null-distributions obtained from shuffled object identities. For minimal variance ratios $Fw,u$⁠, the corresponding values were 18% and 17%, respectively (136 and 130 parcels). To combine the results from both conditions, we used a “prevalence” analysis to determine parcels in which “identity” was represented significantly in a majority of all 16 data sets (prevalence $γ≥0.5$⁠), once again comparing the observed minimal values to the minimal values obtained with shuffled data (red ellipse in Fig. 3C, see Methods for details).

Figure 3A illustrates the 124 parcels identified as significantly “identity-selective” by the prevalence criterion $γ≥0.5$ and Supplementary Figure S3 shows the same information in terms of a sliced brain. Among these were $70$ parcels in the occipital cortex, $29$ in the parietal cortex, $18$ in the fusiform or temporal cortex, and $7$ in the frontal cortex. The average prevalence of identity-selectivity in these parcels was $0.663±0.016$ (mean and S.D.), and the minimal value was $0.58$⁠. As the prevalence criterion (based on 16 data sets) was marginally more conservative than the accuracy criterion (based on 8 data sets), 120 of the 124 parcels were significantly “identity-selective” in terms of both criteria. The four exceptions (identified only by prevalence, but not by accuracy) were Frontal-superior-R 56, Occipital-superior-R 393, Occipital-middle-L 403, and Parietal-superior-R 510. Appendix Table A1 lists the statistical significance of all three criteria for all “identity-selective” parcels.

Overall, there was a pronounced posterior-anterior gradient. Whereas many parcels at the posterior pole of the brain exhibited high classification accuracy, this tended to progressively decrease at more anterior locations (Fig. 3A; Supplementary Fig. S3; Appendix Table A1). To formalize this trend, we assigned 66 of the 124 identity-selective parcels to the 25 topographic visual areas defined by Wang et al. (2015) and, additionally, to the anterior inferior temporal cortex (AIT) and to the inferior frontal cortex (IFC). Supplementary Figure S6 provides an overview of all topographically assigned and non-assigned parcels selective for identity. As illustrated in Figure 8A, this assignment showed that accuracy was comparable in early visual areas (V1-hV4) and in the posterior-ventrolateral regions of the temporal lobe, whereas accuracy was lower in the anterior temporal cortex, the inferior frontal cortex, and in parietal cortical areas.

To assess changes with learning, we separately analyzed five successive and non-overlapping sets of trials (“batches”) with linear classifiers and cross-validation (see Methods for details). Specifically, we established ratios of between- and within-class variance for both object identity (15 classes formed by responses to 15 recurring objects) and for object novelty (2 classes formed by responses to recurring and non-recurring objects, respectively). These two variance ratios measured the neural representation of “identity” and “novelty.”

Variance ratios were converted to z-score values (with respect to the mean and variance of the corresponding shuffle distribution) before being averaged over data sets and/or over parcels. Figure 4A summarizes the results in terms of a grand average over all identity selective parcels. The average identity and novelty ratios were highly significant in all batches (⁠$p<0.001$⁠). Over successive batches, the average identity ratio weakened slightly but significantly (⁠$p<0.05$⁠), whereas the average novelty ratio strengthened considerably, especially between batches $m=1$ and $m=2$ (⁠$p<0.001$⁠).

As expected, it was the between class-variances $SSBidentity$ and $SSBnovelty$ that changed significantly over successive batches $m$ (⁠$p<0.05$ and $p<0.001$⁠, respectively), whereas the within-class variances $SSWidentity$ and $SSWnovelty$ remained essentially the same (⁠$p=n.s.$⁠), as illustrated by Figure 4B. This was owing to the DLDA algorithm, which maintained within-class variance near unity. Nevertheless, over successive batches, the neural representations of recurring objects tended to become slightly more similar to each other, but more dissimilar to the representations of non-recurring objects.

To ascertain that these overall trends hold true also for individual parcels, we carried out more conventional regression analyses of variance ratios $Fm,w,uidentity$ and $Fm,w,unovelty$ over batches $m$⁠, parcels $w$ and data sets $u$⁠. Specifically, we fitted linear mixed-models in order to estimate “rate” parameters $βwidentity$ and $βwnovelty$ for each identity-selective parcel $w$⁠. The results revealed negative rates $βwidentity$ and positive rates $βwidentity$ for almost all parcels, confirming the overall trends in Figure 4C. The variability over parcels was numerically larger for $βwnovelty$ (⁠$0.15±0.1$⁠, mean and S.D.) than for $βwidentity$ (⁠$0.022±0.015$⁠), with both rates weakly correlated (⁠$ρ=0.30$⁠, $p<0.001$⁠). Classification accuracy $αwidentity$ correlated negatively with $βwnovelty$ (⁠$ρ=−0.22$⁠, $p<0.05$⁠) and with $βwidentity$ (⁠$ρ=−​0.74$⁠, $p<0.001$⁠).

To take a closer look at the interaction between “novelty” and “identity,” we divided the identity-selective parcels into “novelty terciles” (high, medium, and low, defined by $βnovelty$⁠) before comparing representations of novelty (⁠$Fnovelty$⁠) and identity (accuracy $α$⁠) (Fig. 5B). The results differed substantially between batches and terciles. In early batches, $Fnovelty$ and $α$ correlated for all terciles, suggesting that initially the representations of non-recurring and recurring objects were linked. However, in successively later batches, this correlation waned in the upper tercile. This may suggest that pronounced representations of non-recurrent objects progressively detached from representations of recurrent objects.