Cell-type-specific microarray data and the Allen atlas: quantitative analysis of brain-wide patterns of correlation and density

For a cell type labeled $t$, we can compute the average correlation with the Allen Atlas in each region (labeled $r$) of the ARA:

For each cell type, the brain regions in the ARA can be ranked by decreasing values of the average correlation $\overline{Corr}(.,t)$ defined in Equation 12. The region with highest average correlation is called the top region by correlation for the cell type labeled $t$. A bar diagram of the average correlations between granule cells (index $t$=20, illustrated in Figure 2(a) of the main text), and the regions of the ’big12’annotation, is shown on Figure 1 (the symbols of the brain regions can be found in Table 1). The top region by correlation for granule cells is the cerebellum.

In the figures derived from the brain-wide correlation profiles (Tables 40–50), the maximal-intensity projections of each correlation profile are supplemented by a section through the top region by correlation.

For a cell type labeled $t$, we can compute the fraction of the total brain-wide density density contributed by voxels of each region labeled $r$ in the ARA:

where Brain Annotation is the set of voxels included in the annotation (for the ’big12’annotation this set consists of the left hemisphere, as can be seen from the projections in Table 1). The brain regions in the ARA are ranked according to the fractions of density defined in Equation 18. For each cell type $t$, the fractional densities supported by the brain regions sum to 1:

The region with highest fraction of density for a given cell type is called the top region by density for this type. A bar diagram of the fractions of the density profile for granule cells (cell type index $t$=20), is shown on Figure 2. The top region by density for granule cells is the cerebellum.

In the figures derived from the brain-wide correlation density profiles (Tables 51–61), the maximal-intensity projections of each density profile are supplemented by a section through the top region by density.

The patchy group of voxels assigned the label ’Basic cell groups and regions’ in the digitized version of the Allen Atlas at a resolution of 200 microns (see Table 1) is not found in the list of brain regions from which the cell-type-specific samples analyzed here were extracted (see Table 2). However, this set of voxels is the top region by correlation and/or density for several cell types, and an inspection of the maximum-intensity projection of the correlation and density profiles for these cell types reveals an anatomical pattern that resembles the brain’s white matter structures, including the arbor vitae. As can be seen in Table 3, most of these cell types were extracted from the cerebral cortex, except the astroglia (sample index $t$=28, [19]), which were extracted from the cerebellum). The sum of the density profiles of these cell types is illustrated on Figure 3. All the cell types whose top region by density is ’Basic cell groups and regions’ are oligodendrocytes, astroglia or astrocytes.

Figures 4 and 5 show results for a class of astrocytes [18] (cell-type index 31) extracted from the cerebral cortex. The brain-wide correlation and density profiles exhibit a pattern resembling white-matter structures, with the most caudal component corresponding to the arbor vitae.

Some of the cell-type-specific samples extracted from the cerebral cortex have a correlation and/or density profile resembling white matter (which is included within the ’Basic cell groups and regions’ in the ARA). See the previous subsection, Table 3 and Figure 3 for a separate analysis of these cell types.

neurons, index 48, and glutamatergic neurons, index 53, both studied in [23]) were extracted from the amygdala (see Tables 64–65). The amygdala in not one of the brain regions in the ’big12’annotation of the ARA. The amygdala is split between the subcortical plate (which is included in the cerebral cortex in the numerical version of the ARA at a resolution of 200 microns) and the olfactory areas [9]. A visual inspection of the correlation and density profiles for both these cell types (Figures 6 and 7) allows to detect a pattern resembling the amygdala, and indeed the cerebral cortex and the olfactory areas rank first and second by the fraction of the density profile they support (see Table 4). It is interesting to examine the fraction of the density profile supported by the various subdivisions of the olfactory areas (according to the ’fine’annotation), especially for glutamatergic neurons, for which olfactory areas support more than 64 % of the total density profile (see Table 5). The Cortical amygdalar area and Piriform-amygdalar area are among the main subregions of the olfactory areas contributing to the density profile of both cell types, which confirms the visual impression of an amygdalar pattern.

Having treated separately the cell types that correlate best or have highest density fits to the white matter (see Table 3), as well as the two amygdalar cell types, we are left with 29 cortical cell types. For 5 of these – all pyramidal neurons – the cerebral cortex is the top region in the ’big12’annotation of the ARA, both by average correlation and by fraction of the density profiles. These pyramidal neurons were all extracted from adult animals, except for the cell type labeled 40 (P14); see Tables 62,63 for age data. See Figure 8–10 for heat maps of the brain-wide correlation and density profiles for these cell types.

Some cell-type specific samples extracted from the cortex as per the anatomical data of Table 2 have the cerebral cortex as top region by correlation, but not as top region by density. There are 8 such cell types (see Table 7), of which 4 have zero density profiles in the left hemisphere. The other four cell types have Olfactory areas or Retrohippocampal region as their top region by density.

Some cortical cell-type-specific samples do not have the cerebral cortex as their top region by correlation or by density. There are 15 such cell types, 7 of which have zero density in the left hemisphere (they are not detected by the linear model, see Table 8). Six cell types have the cortex as their second ranked region by average correlation, while the top region is either the hippocampal region or the retrohippocampal region.

Two cell-type-specific samples that were extracted from the hippocampus (see Table 2). For one of them pyramidal neurons (index 49, studied in [23]), the hippocampal region is the top region both by correlation and by density. We ranked the regions of the ’fine’annotation in the ARA by their contribution to the density profile of this sample (this ranking corresponds to Equation 18, with the region label $r$ running over the $R=94$ regions in the ’fine’annotation). The first two regions are Ammon’s horn (which contributes 48.8 percent of the total density of this cell type), followed by the dentate gyrus (which contributes 25.4 percent of the total density of this cell type). These two regions are the two subregions of the hippocampal region in the ’fine’annotation. Moreover, the data of [23] specify that the sample labeled 49 was taken from Ammon’s horn, which indicates that the brain-wide correlation and density profiles are both consistent with prior biological knowledge for this cell type.

This sample (index 49) is the only one for which the hippocampus is ranked first by density. Another sample (index 3), is ranked first by correlation. This sample was not extracted from the hippocampus, but rather from the cerebral cortex (primary somatosensory area, layer 5). On the other hand, the cerebral cortex is ranked second by its fraction of density defined in Equation 18. This sample has an estimated density profile with values close to zero, except in a few voxels belonging to the olfactory bulb.

The second sample extracted from the hippocampus consists of somatostatin-positive GABAergic interneurons (index 57). The hippocampal region is ranked last by correlations and next-to-last by densities, whereas the first region by average correlation is hypothalamus, and the first region by density is midbrain (see Figure 12 and 14). Visual inspection of the Tables 49,50,60,61 shows that there is a lot of solidarity between the correlation profiles of the cell types labeled with indices between 54 and 64, which are all GABAergic interneurons. For these cell types, with the exception of the cell type labeled 55, whose top region by correlation is the cerebral cortex, the correlations are higher in regions of the brain that are more ventral than the region from which the samples where extracted.

Three cell-type-specific samples (cholinergic neurons – index 13, Drd1+ medium spiny neurons – index 15, and Drd2+ medium spiny neurons – index 16) in the data set were extracted from the striatum.

The pallidum is the region in the ’big12’annotation that supports the highest fraction (43%) of the density of cholinergic neurons (index 13), and the striatum ranks only 7th by correlation and 6th by density (see Table 10 Figures 15 and 18).

On the other hand, medium spiny neurons [19] expressing both dopamine receptor types have striatum as their top region by correlations (see Figures 16, 17, 19, and 20). The striatum is the top region by density for Drd2+ medium spiny neurons, and the second region by density for Drd1+ medium spiny neurons, after the cerebral cortex.

The cell-type-specific sample with index 11 (Cholinergic projection neurons), which is the only sample obtained from the pallidum, has the hypothalamus as its top region by correlation, followed by pallidum (see Figures 21 and 22). Pallidum ventral region, which is the region in the ’fine’ annotation from which the cell-type-specific sample was extracted, ranks 30th out of 94 regions by average correlation (see Table 11). This sugests that cholinergic projection neurons are not very region specific. Moreover, the estimated density of this cell type is zero in the left hemisphere.

Table 12 summarizes results for the cell-type-specific sample extracted from the thalamus (GABAergic Interneurons, PV+, index 59). Thalamus supports less of the density profile than Olfactory areas and Midbrain (it is ranked third by density and fourth by correlation). Moreover, there seems to be little solidarity between the ranking of regions by correlation and by density (see Figures 23 and 24).

Three cell-type-specific samples were extracted from the midbrain. Out of them, two (A9 dopaminergic neurons –index 4– and A10 dopaminergic neurons –index 5) have the midbrain as their top region by density (see Figures 25, 26, 28 and 29). For both of these samples, the midbrain ranks second by correlation, after the hypothalamus.

However, a visual inspection of the projections of correlation and density profiles of Figures 25 and 26 shows some heterogeneity across midbrain, with higher values in the ventral region of it. The ’big12’annotation is too coarse for this heterogeneity to be detected by our ranking procedures. We therefore ranked the regions of the ’fine’annotation are ranked by correlation and density according to Equations 12 and 18. The region that is ranked highest by density for A9 dopaminergic neurons is the ’Substantia nigra, compact part’, which is a subregion of midbrain, and is indeed the finest anatomical label available (see Table 64). This region is ranked second by correlation, the first being ’Ventral tegmental area’ (which is also a subregion of the midbrain).

For A10 dopaminergic neurons the region that is ranked first is the ’Hypothalamus’ (which is the ’generic’ subregion of hypothalamus in the ’fine’annotation, assigned to any voxels that are in hypothalamus but cannot be reliably assigned to a finer subdivision), but the second region is the ’ventral tegmental area, which is a subregion of midbrain, and is indeed the finest anatomical label available from Tables. It is ranked second by correlation, the first being ’Ventral tegmental area’.

Whereas the midbrain supports a majority of the density for A9 dopaminergic neurons and A10 dopaminergic neurons, it supports only 9.8% of the density in Motor Neurons, Midbrain Cholinergic Neurons (whereas pons and medulla support 32.6% and 51.6% respectively, as can be see on Figure 30).

Wheras midbrain supports a majority of the density for A9 dopaminergic neurons and A10 dopaminergic neurons, it supports only 9.8% of the density in Motor Neurons, Midbrain Cholinergic Neurons (whereas medulla and pons support 32.6% and 51.6% respectively). However, it should be noted that very few cell types in our data come from medulla.

One cell-type-specific sample in this study was extracted from the pons (index 51, unpublished data). According to the ’fine’ annotation, it was extracted from the pontine central gray.
Pons is ranked 4th in the ’big12’annotation both by density and correlation (after hypothalamus, midbrain and medulla, which support 28.9%, 26.7% and 12.5% of the density respectively, whereas pons supports 10.4% of the density for this cell type). The pontine central gray is ranked 26th out of 94 regions in the fine annotation. The fraction of the total estimated density of this cell type cumulated by the first 26 regions in the fine annotation is 81.9%.
Conversely, this cell type is the second most important detected in the pons (after the GABAErgic interneurons, PV+, index 64).

The cell-type-specific sample with index 12 (motor neurons, cholinergic interneurons), was extracted from the spinal chord. The closest region in ’big12’is the medulla, which has a refinement into 15 different regions. Medulla is indeed the top region in the ’big12’annotation for this cell type, both by correlation and density (see Figures 33 and 34). The second region in ’big12’by density is the pons (33%), so the medulla and the pons support more than 95% of the estimated density of this cell type. The breakdown of the density amon the regions of the medulla is shown in Table 16.

Out of the 11 samples that were drawn from the cerebellum, 7 have the cerebellum as their top region by correlation, and 4 have the cerebellum as their top region by density (all of which also have the cerebellum as their top region by correlation, see 17). All the samples that have the cerebellum as their top region were indeed taken from the cerebellum. See Table 3 for the Mature oligodendrocytes (index 28), whose density profiles follows a white-matter pattern that includes the arbor vitae. See Figure 35 for a class of Purkinje cells [23] and Figure 36 for a class of mature oligodendrocytes [19], both extracted from the cerebellum. Their correlation and estimated density patterns are indeed mostly localized in cerebellum.

A remarkable class of Purkinje cells does not have the cerebellum as its top region by density. See Figure 37 for a class of Purkinje cells (index 52, unpublished) that correlates best with the Allen Atlas both in the thalamus and in the cerebellum, but that fits based on density coefficients only in the thalamus. This indicates that thalamus must contain cell types whose gene expression profile is closest to Purkinje cells in the present data set, but that thalamus is not sampled in enough detail by our microarray data set for these cell types to be distinguished from this class of Purkinje cells.

We looked for families of cell-type-specific samples that are very similar to each other. We computed the pairwise type-by-type correlation matrix between centered cell-type specific transcription profiles:

The closer the entry $typeCorr(t,t^{\prime })$ is to 1, the more similar the two cell types $t$ and $t^{\prime }$ are. We applied the thresholding procedure that was used for co-expression matrices of genes in the brain in [31]. For each value of $\tau$ in the interval $[-1,1]$, a threshold can be applied to the matrix $typeCorr$ by putting to zero all entries lower than $\tau$:

Applying Tarjan’s algorithm [Tarjan] to the thresholded matrix ${typeCorr}_{\rho }$ produces a partition of the $T$ cell types into strongly connected components at this value of the threshold. At $\tau =1$, all the cell types are disconnected (unless there are cell types with exactly the same value for all genes, which is not the case in our data set). At $\tau =-1$, there is only one connected component, and all cell types are connected. When $\tau$ decreases from $1$ to $-1$, the more similar cell types group into connected components. We are interested in pairs of cell types that are connected at high values of $\tau$ and do no not connect with other cell types for the lowest possible value of $\tau$, when $\tau$ is decreased (such cell types labeled $t$ and $t^{\prime }$ are the closest pairs of cell types to the degenerate situation described by Equation 4.1.1.

Consider a pair of indices $(t_1,t_2)$ singled out by the above analysis. Let us introduce the notation ${\hat{C}}^i$ for the type-by-gene matrix obtained from $C$ by leaving the $i$-th row out, we computed the corresponding density profiles denoted by ${\hat{\rho }}^{t_1}$ and ${\hat{\rho }}^{t_2}$:

The pairs of cell types singled out by the above-described analysis are the following:
$\left(i\right)$ indices $(15,16)$, which are both medium spiny neurons (Drd1+ and Drd2+ respectively),
$\left(ii\right)$ indices $(4,5)$, which are both dopaminergic neurons (A9 and A10 respectively),
$\left(iii\right)$ indices $(2,3)$, which are both pyramidal neurons (but are not detected by the model except at a few voxels in the olfactory bulb).

• The pair of medium spiny neurons (labeled by indices $t_1=15$ and $t_2=16$). Having computed the refitted profiles ${\hat{\rho }}^{15}$ and ${\hat{\rho }}^{16}$, which contain one medium spiny neurons instead of two, we can plot the density profiles oh the remaining medium spiny neurons in each of the refitted models (${\hat{\rho }}_{16}^{15}$ and ${\hat{\rho }}_{15}^{16}$, Figures 38(b) and 38(c) respectively), as well as the sum of the densities of the two samples of medium spiny neurons in the original model (i.e. ${\rho }_{15}+{\rho }_{16}$, Figure 38(c)). The three plots are hard to distinguish from each other by eye (see Figure 38) and present the same characteristic concentration of signal in the striatum. This supports the conjectures 24 and 24, as each of the medium spiny neurons seems to inherit the signal of the pair labeled by indices $15$ and $16$ when the model is refitted, without picking up any significant other signal.
  

  

  Computing the correlation matrices ${\hat{\Gamma }}^{15}$ and ${\hat{\Gamma }}^{16}$ and reading off the entries

  $${\hat{\Gamma }}_{16}^{15}=0.9997,{\hat{\Gamma }}_{15}^{16}=0.9992$$

  (35)

  confirms the visual impression of Figure 38, due to definitions 30 and 33.
  

  Moreover, plotting the matrices ${\hat{\Gamma }}^{15}$ and ${\hat{\Gamma }}^{16}$ as heat maps (Figures LABEL:projectionsPairfit_15_16b and Figures LABEL:projectionsPairfit_15_16c) shows that the entries are roughly symmetric, and that the two matrices resemble the type-by-type matrix of $\Gamma$ of pairwise correlations in the original model (Figure LABEL:projectionsPairfit_15_16a ), apart from the “cross” of zeroes at the omitted indices (row 15 and column 15 on Figure LABEL:projectionsPairfit_15_16b, row 16 and column 16 on Figure LABEL:projectionsPairfit_15_16c). Moreover, the diagonal coefficients of ${\hat{\Gamma }}^{15}$ and ${\hat{\Gamma }}^{16}$ are visibly close to 1 (their sorted values are plotted on Figure 39), with the exception of the cholinergic projection neurons (labeled $t=11$), whose density profile ${\widehat{{\rho }^{15}}}_{11}$ has only a 0.133 correlation coefficient with ${\rho }_{11}$. However, this cell type has very little signal in any of the models (${\rho }_{11}$ has only 7 voxels with positive density, ${\hat{\rho }}_{11}^{15}$ has 8 and ${\hat{\rho }}_{11}^{16}$ has 7), and these three profiles all represent less than ${10}^{-8}$ times the sum of all the densities in the respective models:

  $$\frac{{\sum }_v{\rho }_{11}\left(v\right)}{{\sum }_t{\sum }_v{\rho }_{11}\left(v\right)}=1.53\times {10}^{-9},\frac{{\sum }_v{\hat{\rho }}_{11}^{15}\left(v\right)}{{\sum }_t{\sum }_v{\hat{\rho }}_{11}^{15}\left(v\right)}=6.65\times {10}^{-9},\frac{{\sum }_v{\hat{\rho }}_{11}^{16}\left(v\right)}{{\sum }_t{\sum }_v{\hat{\rho }}_{11}^{16}\left(v\right)}=1.53\times {10}^{-9}.$$

  (36)

  The low value of the correlation for cell-type index 11 is therefore compatible with the claim that density profiles are generally stable when replacing the pair of medium spiny neurons by just one of them, and that the remaining medium spiny neuron in the fitting panel inherits the sum of the density profiles predicted in the original model.

  \@float

  figure\end@float

  

• The pair (A9 dopaminergic neurons, A10 dopaminergic neurons) labeled by indices $t_1=4$ and $t_2=5$. We computed the refitted profiles ${\hat{\rho }}^4$ and ${\hat{\rho }}^5$, which contain one sample of dopaminergic neurons instead of two. Again we can plot the density profiles of the remaining dopaminergic neuron in each of the refitted models (${\hat{\rho }}_5^4$ and ${\hat{\rho }}_4^5$, Figures 40(b) and 40(c) respectively), as well as the sum of the densities of A9 dopaminergic neurons and A10 dopaminergic neurons in the original model (i.e. ${\rho }_4+{\rho }_5$, Figure 40(c)). The three plots show the highest fraction of their signal in the midbrain (hence the sections through mibrain in all three rows of Figure 40). generally speaking they look quite similar. Again this supports the conjectures 24 and 24, in the case of $t_1=4$ and $t=5$.
  

  

  Computing the correlation matrices ${\hat{\Gamma }}^4$ and ${\hat{\Gamma }}^5$ and reading off the entries The correlations conjectured in Equations 30 and 33 are less close to 1 than in the case of medium spiny neurons, but still

  $${\hat{\Gamma }}_5^4=0.87,{\hat{\Gamma }}_4^5=0.96,$$

  (37)

  and one can check that in each of the two refitted models, the sum of profiles ${\rho }_4+{\rho }_5$ has higher correlation with the remaining dopaminergic neurons than any cell type in the original fitting panel:

  $${\hat{\Gamma }}^{t_1}(t_2,t_2)={max}_{1\le t\le T}{\hat{\Gamma }}^{t_1}(t_2,t),{\hat{\Gamma }}^{t_2}(t_1,t_1)={max}_{1\le t\le T}{\hat{\Gamma }}^{t_2}(t_1,t).$$

  (38)

  As in the case of medium spiny neurons, plotting the matrices ${\hat{\Gamma }}^4$ and ${\hat{\Gamma }}^5$ as heat maps (Figures LABEL:projectionsPairfit_4_5b and Figures LABEL:projectionsPairfit_4_5c) shows that the entries are roughly symmetric, and that the two matrices resemble the type-by-type matrix of $\Gamma$ of pairwise correlations in the original model (Figure LABEL:projectionsPairfit_4_5a ), apart from the “crosses” of zeroes (at row 4 and column 15 on Figure LABEL:projectionsPairfit_4_5b, row 16 and column 16 on Figure LABEL:projectionsPairfit_4_5c). Moreover, the diagonal coefficients of ${\hat{\Gamma }}^4$ and ${\hat{\Gamma }}^5$ are all close to 1 (their sorted values are plotted on Figure 41). Again we conclude that refitting the model to a fitting panel containing one dopaminergic neuron instead of two leads to stable density profile, with the remaining dopaminergic neuron inheriting most of the sum of the two original density profiles ${\rho }_4$ and ${\rho }_5$.

  \@float

  figure\end@float

  

• The pair of pyramidal neurons labeled by indices $t_1=3$ and $t_2=3$. We computed the refitted densities ${\hat{\rho }}^2$ and ${\hat{\rho }}^3$ and repeated the above analysis. Figures LABEL:projectionsPairfit_2_3, 43 and 42 confirm conjectures 24 and 24. However, they are included mostly for completeness, as the estimated densities of these two pyramidal neurons are very low in the original model (only 52 voxels , about 1 in 1000, have a positive density in the sum ${\rho }_2+{\rho }_3$). They stay low in the refitted model (${\hat{\rho }}_3^2$ has only 49 voxels with non-zero density and ${\hat{\rho }}_2^3$ only 44). The high diagonal correlation coefficients of Figure 43 show that the rest of the estimated densities hardly change when refitting.

  

  \@float

  figure\end@float