GWAS-based pathway analysis differentiates between fluid and crystallized intelligence.

Cognitive abilities vary among people. About 40-50% of this variability is due to general intelligence (g), which reflects the positive correlation among individuals' scores on diverse cognitive ability tests. g is positively correlated with many life outcomes, such as education, occupational status and health, motivating the investigation of its underlying biology. In psychometric research, a distinction is made between general fluid intelligence (gF) - the ability to reason in novel situations - and general crystallized intelligence (gC) - the ability to apply acquired knowledge. This distinction is supported by developmental and cognitive neuroscience studies. Classical epidemiological studies and recent genome-wide association studies (GWASs) have established that these cognitive traits have a large genetic component. However, no robust genetic associations have been published thus far due largely to the known polygenic nature of these traits and insufficient sample sizes. Here, using two GWAS datasets, in which the polygenicity of gF and gC traits was previously confirmed, a gene- and pathway-based approach was undertaken with the aim of characterizing and differentiating their genetic architecture. Pathway analysis, using genes selected on the basis of relaxed criteria, revealed notable differences between these two traits. gF appeared to be characterized by genes affecting the quantity and quality of neurons and therefore neuronal efficiency, whereas long-term depression (LTD) seemed to underlie gC. Thus, this study supports the gF-gC distinction at the genetic level and identifies functional annotations and pathways worthy of further investigation.

A large proportion (40–50%) of interindividual variability in cognitive abilities is due to general intelligence (g), a quantitative trait that reflects the positive correlation among an individual's scores on diverse cognitive ability tests (Deary 2012; Spearman 1904). A high g score is associated with many favorable life outcomes (Deary & Batty 2011). Psychometric research distinguishes between general fluid intelligence (gF) – the ability to reason in novel situations – and general crystallized intelligence (gC) – the ability to apply acquired knowledge and learned skills (Carroll 1993; Cattell 1963). Although gF and gC are correlated at least 50% based on twin studies (Wainwright ) and more so early and late in life (Li ), developmental and cognitive neuroscience studies largely support the distinction between them. In normal aging, gF declines earlier and more rapidly than gC (Craik & Bialystok 2006; Salthouse 2004); in development, measures of verbal (gC) and non-verbal (gF) intelligence correlate differently with changes in brain structure (Ramsden ). Furthermore, gF is closely associated with fronto-parietal and anterior cingulate brain networks, whereas gC is dependent on posterior frontal and temporal regions (Barbey ; Glascher ; Jung & Haier 2007; Woolgar ). Finally, at the population level, large gains in performance have been observed for tests that are strongly associated with gF, but not with gC (Flynn 2007).

More than half of the variability in intelligence tests is attributable to additive genetic effects (Deary ; Lee ; Plomin & Spinath 2004). In a recent genome-wide association study (GWAS), 40% and 51% of the phenotypic variability in gC and gF, respectively, could be accounted for by genetic variants in linkage disequilibrium (LD) with common single nucleotide polymorphisms (SNPs) (Davies ), providing a lower-bound estimate of the narrow-sense heritability of these traits. Furthermore, using only SNP data, ∼1% of the variance in intelligence test scores in a sample could be predicted in an independent sample (Davies ). However, at the single-marker level, no robust genetic association with intelligence has yet been published, consistent with the observation that the effect of individual SNPs may be too weak to be detected in complex polygenic traits using the classical genome-wide P ≤ 5 × 10−8 threshold (Wellcome Trust Case Control Consortium 2007) on the sample sizes currently available (International Schizophrenia Consortium ). Therefore, relaxing the significance criteria and exploiting the polygenic signal by going beyond the traditional single-marker approach to gene- and pathway-based methods may offer more power (Neale & Sham 2004) and insight into the biological processes underlying these traits (Wang ).

The polygenic architecture of gF and gC was previously confirmed using the Cognitive Aging Genetics in England and Scotland (CAGES) cohort as the discovery sample and the Norwegian Cognitive NeuroGenetics (NCNG) adult lifespan sample for replication (Davies ). Here, using the same datasets, we combine single-marker, gene- and pathway-based approaches to characterize the genetic architecture of gF and gC with respect to known biological processes.

Materials and methods

GWASs: genotypes and phenotypes

The genotype and phenotype protocols and the samples have been described previously (Davies ; Espeseth ). All participants gave written consent before the study started. All procedures were conducted according to the tenets of the Declaration of Helsinki and approved by the relevant Research Ethics Committees. An overview of the genotype and phenotype protocols in the discovery and replication samples is provided below, with further details available in the original publications.

Discovery GWAS

The final NCNG GWAS consisted of 554 225 SNPs genotyped in a homogenous Norwegian sample of 670 individuals (457 females), ranging from 18 to 79 years of age (M = 47.6; SD = 18.3) (Espeseth ). Participants completed a battery of psychometric tests, assessing general cognitive ability, memory, attention and speed of processing. The protocol was approved by the Regional Committee for Medical and Health Research Ethics, Southern Norway (project ID: S-03116). DNAs were genotyped on the Illumina Human610-Quad Beadchip. Quality control was performed with the ‘check.marker’ function of the R package GenABEL (Aulchenko ). Individuals were excluded based on relatedness (‘ibs.threshold’ = 0.85), heterogeneity, unresolved sex discrepancies and call rate ≤ 0.97. Population structure was assessed by multidimensional scaling (MDS) analysis (using 100 000 random SNPs), removing outlying samples with possible recent non-Norwegian ancestry. No additional adjustment for population structure was performed because of the homogeneity of the sample (Espeseth ). Single nucleotide polymorphisms with a call rate ≤ 0.95, minor allele frequency ≤ 0.01 and Hardy–Weinberg equilibrium (exact test) P-value ≤ 0.001 were excluded.

Replication GWAS

The CAGES sample consists of five cohorts, the Lothian Birth Cohorts of 1921 (Deary , 2009b) and 1936 (Deary ) (LBC1921, LBC1936), the Aberdeen Birth Cohort of 1936 (ABC1936) (Deary , 2009b) and the Manchester and Newcastle Longitudinal Studies of Cognitive Aging Cohorts (Rabbitt ). The final GWAS dataset consisted of 549 692 SNPs genotyped in 3511 healthy individuals (2115 females) with detailed cognitive ability measurements taken in middle to older adulthood (mean age ranged from 64.6 in the ABC1936 to 79.1 in LBC1921; overall age range: 44–93) (Davies ). Ethical approval was obtained from the relevant Research Ethics Committees. Individuals were excluded based on unresolved gender discrepancy, coefficient of relatedness > 0.025 (based on 549 692 autosomal SNPs), call rate ≤ 0.95 and evidence of non-Caucasian descent as determined by MDS analysis (Davies ). Single nucleotide polymorphisms were included if they met the following conditions: call rate ≥ 0.98, minor allele frequency ≥ 0.01 and Hardy–Weinberg equilibrium test with P ≥ 0.001. Population structure was assessed by MDS analysis, and four MDS components were fitted as covariates to correct for any population stratification that might be present.

Cognitive phenotypes

Different measures of gF and gC were employed for each of the five CAGES samples and the NCNG sample (Davies ). In general, the different tests of gF aimed to assess each individual's capacity to reason logically and solve problems in novel situations, relatively independently of acquired knowledge. gF employs aspects of processing speed, attention, memory and executive function. The gC tests were vocabulary-based, assessing each individual's semantic knowledge. All measures were corrected for sex (with the exception of the CAGES Manchester and Newcastle gF, which was derived separately for males and females) and age. The standardized residuals were then extracted and used as the trait measures in all subsequent analyses. In the NCNG, of the 670 individuals, 629 and 643 had scores for gF and gC, respectively. In the CAGES, of the 3511 individuals, 3400 and 3482 had scores for gF and gC, respectively.

Analyses

Construction of gene lists

A ‘bottom-up’ approach (Liu ) was undertaken. Both single-marker and gene-based association analyses were performed, testing SNPs and genes individually first in order to construct a list of genes that would then be subjected to pathway analysis, via Ingenuity Pathway Analysis (IPA), to identify over-represented functions and/or pathways. For each of gF and gC, genetic factors that showed evidence of association in the NCNG were identified and filtered further based on evidence of replication in the CAGES. Single-marker and gene-based statistical methods and thresholds were selected and applied, as appropriate, with the aim of controlling the Type II (false-negative) as well as the inversely related Type I (false-positive) error rates, incorporating instead prior biological knowledge to the interpretation of the findings and generating testable hypotheses for further investigation (Lieberman & Cunningham 2009; Williams & Haines 2011).

The LDsnpR tool (Christoforou ) was used to (1) annotate individual SNPs and (2) assign SNPs to genes for the gene-based analyses. Single nucleotide polymorphisms were assigned to genes (Ensembl 54 definitions) if they were located physically within the boundaries of the gene (±10 kb) or if they were in high LD [r2 ≥ 0.80 based on HapMap CEU (http://hapmap.ncbi.nlm.nih.gov/)] with another SNP located within the boundaries of the gene (±10 kb).

Single-marker analyses

In the NCNG sample, association between individual SNP alleles and cognitive phenotypes was tested using linear regression analysis (1 d.f. coefficient t-test), as implemented in PLINK (Purcell ). The CAGES GWAS data comprised the inverse variance weighted model ‘meta P-values’ produced in the original CAGES five-sample meta-analysis (i.e. as reported by Davies ). A meta-analysis of the CAGES and NCNG single-marker allele P-values was performed on the overlapping SNPs, using an inverse variance weighted model, as implemented in METAL (Willer ).

All SNPs with an asymptotic nominal P-value of ≤0.05 in the NCNG were then mined for replication in the CAGES GWAS data. As advised in Konig (2011), positive single-marker replication was determined on the basis of P ≤ 0.05 and the same direction of effect in both the NCNG and CAGES samples, resulting in a meta-analysis P-value that is more significant than either of the original P-values. Using LDsnpR, all SNPs that replicated were annotated, where possible, with a gene name or Ensembl identifier and taken forward for pathway-based analysis by IPA.

Gene-based analyses

All 554 225 SNPs in the NCNG GWAS were assigned, where possible, to genes using LDsnpR. For each gene containing an SNP, association was assessed by applying PLINK's permutation-based set test on the LDsnpR-generated sets. Gene-based statistics generated by permutation-based methods automatically account for potential confounding factors, such as LD structure and gene length or SNP number, through the generation of an empirical null distribution (Liu ). As the true underlying genetic architecture of these traits is unknown, three different models were tested to obtain a gene-based association score (Lehne ):

The ‘minimum P-value’ (minP) model, which assigns to each gene the association statistic, or P-value, of the most significant SNP in the gene, assumes that a single SNP within the gene contributes to the phenotype.

The ‘all P-values’ model, which assigns to each gene the mean association statistic of all the SNPs in the gene, assumes that all or most SNPs within the gene contribute to the phenotype.

The ‘P ≤ 0.05 threshold’ model, which assigns to each gene the mean association statistic of all SNPs within the gene that have a P ≤ 0.05, assumes that only a few or a subset of SNP contribute to the phenotype.

The P-values were computed based on 10 000 permutations, and the lowest of the three permuted P-values was retained as the gene-based P-value for each gene. No additional correction was applied to the gene-based score to account for testing the three models because of high correlation between the tests. Genes with gene-based permuted P ≤ 0.05 were mined for replication in the CAGES.

In order to preserve the Type II error rate, the replication criteria in the CAGES were relaxed to include any gene that contained at least one SNP with a ‘meta-P’ ≤ 0.05, ensuring that significance under any genetic architecture was captured. In the NCNG, for example, all genes with a minimum gene-based permuted P ≤ 0.05 contained at least one SNP with a P ≤ 0.05 (data not shown). This replication criterion of at least one SNP with a ‘meta-P’ ≤ 0.05 is equivalent to the unadjusted minimum P-value approach. Thus, the CAGES SNPs were first assigned to genes using LDsnpR, as described above, and the genes were scored using the minimum P-value approach, assigning to the gene the P-value of the most significant SNP in that gene (without further adjustment). All genes that were selected from the NCNG and that also replicated in the CAGES (i.e. meta-P ≤ 0.05) were taken forward for analysis with IPA. Finally, after the relevant gene lists were pruned for LD on the basis of their association signals and position, the hypergeometric distribution, as implemented in www.geneprof.org/GeneProf/tools/hypergeometric.jsp, was used to assess whether the number of genes replicating in the CAGES was greater than expected by chance.

Construction of gene lists – stringent replication criteria

As the aforementioned single-marker and gene-based replication criteria in the CAGES are prone to inflating the Type I error rate due to gene-length bias, the robustness of the subsequent main IPA findings was assessed by re-running IPA on gene lists constructed using more stringent statistical thresholds and methods. First, for the single-marker analysis, only genes tagged by the SNPs that met the previously described replication criteria and resulted in a meta-analysis P ≤ 1 × 10−5, as supported by the National Human Genome Research Institute GWAS catalog (Welter ), were taken forward for IPA analysis. For the gene-based analysis, replication in the CAGES was determined using two gene-scoring methods, which accounted for the number of SNPs per gene and, critically, for LD without requiring genotype-level data. The first method scored each gene with the most significant P-value (i.e. the minP approach), adjusting for the number of SNPs using a modified Sidak's correction (Saccone ). This approach has been shown to perform as well as regression-based methods (Christoforou ; Segre ). The second approach was comparable to the regression-based ‘all P-values’ model described above, combining all SNPs mapped to the gene using Brown's approximation method, as implemented in PLINK (Moskvina ; Purcell ). The NCNG genotype data was used to estimate the LD between SNPs. Thus, the genes with gene-based permuted P ≤ 0.05 in the NCNG, which met a P-value of ≤ 0.05 with either of these two methods in the CAGES, were subjected to IPA analysis.

Ingenuity Pathway Analysis

Genes that showed evidence of association in the NCNG and of replication in the CAGES based on either the single-marker and/or the gene-based analysis were selected for pathway analysis with IPA (Ingenuity® Systems, www.ingenuity.com) to identify the most relevant biological functions and pathways. Ingenuity Pathway Analysis gathers information from published data from ∼3600 peer-reviewed journals regarding relationships between genes and proteins. The information is updated weekly, manually curated and stored within the Ingenuity® Knowledge Base, which is then queried during an analysis to identify specific biological functions, or ‘function annotations’, and pathways enriched within the submitted gene list.

In this study, the Ensembl 54 gene identifiers were uploaded into IPA and mapped, if possible, to their corresponding object in the Ingenuity® Knowledge Base (Genes Only). A ‘Core Analysis’ was performed, including both direct and indirect relationships and using all available data sources in all species available. Molecules and/or relationships that were either experimentally observed or predicted with high confidence were considered. For the main analyses, which aimed to differentiate gF and gC as it relates, most relevantly, to the nervous system, tissues and primary cells were restricted to those of the nervous system, and cell lines were restricted to those of the central nervous system. Additional analyses were performed to evaluate the robustness and relevance of the main findings, including (1) re-running the IPA having excluded genes known to be extensively studied to assess the possibility of publication bias driving the main findings and (2) re-running the IPA having included all tissues and cell lines, except those specific to cancer, to assess the relevance of the main findings in the context of more global annotation.

In particular, we focused on the significant function annotations. The function annotations are organized based on Ingenuity® Ontology, which consists of a manually built and maintained hierarchical data structure of hundreds of thousands of unique classes. As part of the ‘Core Analysis’, IPA performs a Fisher's exact test (FET) to determine whether the submitted gene list consists of genes of a particular function annotation more than expected by chance, given the proportion of genes of that particular function annotation in the entire Ingenuity Knowledge Base. It also uses the Benjamini–Hochberg (BH) multiple-testing correction method to adjust for the number of ontologies tested, providing a false discovery rate for a particular function annotation. In this study, all functional annotations with an FET P ≤ 0.05 (and corresponding BH P ≤ 0.25) were considered significant.

Results

Construction of gene lists

Single-marker analyses

A total of 554 225 SNPs were tested for allelic association to gF (N = 629 individuals) and gC (N = 643 individuals) in the NCNG sample. As these GWAS data results were not presented in the previous study (Davies ), the traditional SNP-level diagnostics and results are provided in the Supporting Information (Fig. S1 and Dataset S1). The genomic inflation factor was 1 for both gF and gC (Fig. S1), indicating that the data did not suffer from population stratification or other systematic bias. When mining the CAGES data, 816 SNPs in gF and 884 SNPs in gC met the replication criteria of P ≤ 0.05 and the same direction of effect (Dataset S1). Using the LD-based binning approach implemented in LDsnpR (Christoforou ), 481 (59%) of the gF SNPs and 549 (62%) of the gC SNPs were assigned to Ensembl 54 genes. 159 gF SNPs and 160 gC SNPs mapped to more than one gene on the basis of LD. As it is not possible to identify the true source of the association on the basis of the SNP P-value alone (Christoforou ,2012b), all genes were retained for pathway analysis. The significant function annotations or pathways were subsequently manually examined to ensure that they did not consist of clusters of genes representing the same genetic association signal. Thus, 503 and 530 Ensembl 54 genes for gF and gC, respectively (Dataset S1), were taken forward for IPA analysis. A total of 171 (34%) of the gF genes and 212 (40%) of the gC genes were implicated by more than one replicated SNP.

Gene-based analyses

A gene-based approach was also used to identify candidate genes for pathway analysis, allowing for locus heterogeneity and the aggregation of multiple weaker association signals. In the NCNG, of the 34 109 eligible Ensembl 54 gene entries, 2698 and 2615 met the nominal significance threshold of set-based permuted P ≤ 0.05 in gF and gC, respectively (Dataset S2). Of these, 841 gF genes and 920 gC genes had a minimum P-value gene-based score of ≤ 0.05 in the CAGES (Dataset S2). After the relevant gene lists were first pruned for LD based on position and association signal, resulting in 607 and 652 relatively independent association signals for gF and gC, respectively, the number of genes showing evidence of replication was more than expected by chance for gC (hypergeometric P = 4.5 × 10−4), but not for gF (hypergeometric P = 0.25). The 841 gF genes and 920 gC genes were added to the list of genes identified by single-marker analysis and taken forward for IPA.

Figure 1 summarizes the number of SNPs and genes identified in each analysis and the relevant corresponding overlaps. A total of 1182 genes for gF and 1294 genes for gC were identified through both the single-marker and gene-based analyses and subjected to IPA analysis. Of these genes, 853 for gF and 893 for gC were identified by their Ensembl Gene identifier in the IPA database (i.e. ‘IPA Ready’ genes; Dataset S3) and were thus available for pathway analysis. A total of 117 genes were in common between the two traits.

Figure 1

Summary of SNP- and gene-based analyses and corresponding numbers of genes identified for IPA for both gF (left) and gC (right). The boxed numbers in the centre represent the number of genes identified as a result of the SNP-based analyses (from top, down), the gene-based analyses (from bottom, up), the resulting total number of unique genes submitted to IPA and the number of relevant overlapping genes (in parentheses). * indicates that the GWAS SNPs were assigned to genes by LDsnpR.

IPA's ‘Core Analysis’ was run on the two lists of genes, restricting the analysis to tissues and cell lines related to the nervous system. Function annotations that were significantly enriched with an FET P ≤ 0.05 and a BH P ≤ 0.25, which account for all of the function annotations tested, were declared significant and evaluated. For both traits, function annotations involved in the biological function ‘Nervous System Development and Function’ were the most common, accounting for >25% of the function annotations. These were followed by annotations involved in ‘Cell-to-Cell signaling’ and ‘Cellular Assembly and Organization’ in both gF and gC (∼8%) (Dataset S3). However, it is important to note that some function annotations are categorized under multiple biological functions, resulting in redundancy in the data (e.g. ‘axonogenesis’ appears under ‘Nervous System Development and Function’, ‘Cell Morphology’, ‘Cellular Function and Maintenance’, ‘Cellular Assembly and Organization’ and ‘Tissue Development’). When considering only the non-redundant function annotations, specifically those identified on the basis of at least two genes from the inputted gene list, 85 and 54 function annotations were identified as significantly enriched in gF and gC, respectively (Dataset S3).

Only four function annotations were enriched in both gF and gC, including ‘microtubule dynamics’ and ‘chemotaxis of neurons’, which could be generally categorized as architectural features of neuronal development. The set of genes leading to the enrichment of these overlapping annotations in gF was different to that in gC (Table1), except for the ‘migration of GABAergic neurons’ annotation which was attributed to the same two genes in gF and gC.

Table 1

Significant IPA function annotations (FAs) in both gF and gC

Function	FA	gF			gC
		FET P-value	BH P-value	Genes	FET P-value	BH P-value	Genes
Chemotaxis	Chemotaxis of neurons	0.025	0.18	BDNF, GFRA1, RGS3	0.0037	0.17	EPHB2, GFRA1, SEMA3A, SLIT2
Formation	Formation of dendrites	0.028	0.18	ACHE, BCL11A, BDNF, CTNNA2, GRIN3A, NRG1, PRKG1, SGK1, SYNE1	0.035	0.22	DAB1, DSCAM, KLF7, KNDC1, MAP2, NRG1, PPP1R9B, RAC2, SEMA3A
Organization	Microtubule dynamics	0.030	0.18	ACHE, ATXN10, BAX, BCL11A, BCL2, BDNF, BSN, CDH1, CNTN4, CNTNAP2, CTNNA2, CTNND2, CYP19A1, DISC1, EGFR, GDA, GRIN3A, GSN, IGF1R, KATNB1, MAPK8, NFIB, NRG1, PRKG1, PTPRM, RND1, SGK1, SLC18A3, SLIT1, SLIT3, SNCA, SYNE1, UHMK1	0.032	0.22	ATXN2, BBS10, CDH4, CNTN4, DAB1, DCC, DSCAM, EPHB1, EPHB2, GPM6A, KLF7, KNDC1, LAMB1, LRRC7, MAP2, MBP, MTOR, NRG1, PCDH15, PLD1, POU3F1, PPP1R9B, PRICKLE2, PRKCA, PTPRM, RAC2, RIT2, RUFY3, SEMA3A, SLIT1, SLIT2, TNFRSF12A, TRPV4, VANG
Migration	Migration of GABAergic neurons	0.012	0.18	GFRA1, SLIT1	0.013	0.22	GFRA1, SLIT1

For each general function, the specific FA is shown together with the corresponding FET P-value, the BH multiple-testing corrected P-value and the genes responsible for the enrichment signal in the FA. Individual genes that contribute to the enrichment of a particular FA in both gF and gC are in bold.

The most statistically significant functional annotations were identified in gC. ‘Synaptic depression’ ranked at the top (FET P = 2.9 × 10−6; BH P = 0.0015), with the related ‘long-term depression’ in general (LTD; FET P = 2.0 × 10−5; BH P = 0.0052) and LTD of specific cells (Table2), all of which are highly related as indicated by the composite genes. ‘Guidance of axons’ and ‘schizophrenia’ were also significantly enriched in gC (FET P = 4.8 × 10−4 and 0.042, respectively), with the former also ranking third (Table2 and Dataset S3).

Table 2

Top 25 IPA FAs for gC

Function	FA	FET P-value	BH P-value	No. of genes	Genes
Synaptic depression	Synaptic depression	2.9 × 10−6	0.0015	13	ADCY5, ADCY8, CNR1, DRD5, EPHB1, EPHB2, GRM7, LRRC7, MTOR, PRKCA, RYR3, ST8SIA4, SYNJ1, (ITSN1)
LTD	LTD	2.0 × 10−5	0.0052	11	ADCY5, ADCY8, CNR1, DRD5, EPHB2, GRM7, LRRC7, MTOR, PRKCA, RYR3, ST8SIA4
Guidance	Guidance of axons	4.8 × 10−4	0.061	16	ALCAM, ANK3, CDH4, CNTN4, DCC, EPHB1, EPHB2, ERBB4, EXT1, GLI3, KLF7, NFASC, PLXNA2, SEMA3A, SLIT1, SLIT2
LTD	LTD of cells	4.8 × 10−4	0.061	6	ADCY5, CNR1, DRD5, PRKCA, RYR3, ST8SIA4
LTD	LTD of brain cells	7.2 × 10−4	0.075	5	ADCY5, CNR1, PRKCA, RYR3, ST8SIA4
LTD	LTD of neurons	0.0013	0.11	5	ADCY5, CNR1, DRD5, RYR3, ST8SIA4
Assembly	Assembly of olfactory cilia	0.0020	0.13	3	BBS10, PCDH15, VANGL2
LTD	LTD of corticostriatal neurons	0.0023	0.13	2	ADCY5, CNR1
Positioning	Positioning of cholinergic neurons	0.0023	0.13	2	SLIT1, SLIT2
Contact repulsion	Contact repulsion	0.0034	0.17	3	DCC, SEMA3A, SLIT2
Chemotaxis	Chemotaxis of neurons	0.0037	0.17	4	EPHB2, GFRA1, SEMA3A, SLIT2
Development	Development of optic chiasm	0.0068	0.22	2	SLIT1, SLIT2
Size	Size of lateral cerebral ventricle	0.0068	0.22	2	ANK2, NRG1 (includes EG:112400)
Apoptosis	Apoptosis of sensory neurons	0.0076	0.22	3	CDKN2D, HIPK2, KLF7
Metabolism	Fatty acid metabolism	0.010	0.22	5	DAB1, ERBB4, NRG1, SEMA3A, ST8SIA1
Abnormal morphology	Abnormal morphology of nervous system	0.013	0.22	2	CNR1, UHRF1
Formation	Formation of oligodendrocytes	0.013	0.22	2	ERBB4, NRG1
Guidance	Guidance of thalamocortical axons	0.013	0.22	2	SLIT1, SLIT2
Migration	Migration of GABAergic neurons	0.013	0.22	2	GFRA1, SLIT1
Pathfinding	Pathfinding of axons	0.013	0.22	2	DCC, EXT1
Chemotaxis	Chemotaxis of cells	0.014	0.22	5	EPHB2, FPR2, GFRA1, SEMA3A, SLIT2
Development	Development of diencephalon	0.014	0.22	3	SIM2, SLIT1, SLIT2
LTD	LTD of synapse	0.018	0.22	3	CNR1, GRM7, MTOR
Synthesis	Synthesis of fatty acid	0.020	0.22	4	DAB1, ERBB4, NRG1, SEMA3A
Synthesis	Synthesis of lipid	0.020	0.22	6	CNR1, DAB1, ERBB4, NRG1, SEMA3A, ST8SIA1

For each general function, the specific FA is shown together with its corresponding FET P-value, the BH multiple-testing corrected P-value, the number of genes and the symbols of the genes responsible for the enrichment signal in the FA. The functions and FAs that remained significant in the IPA analysis of the genes that passed the more stringent criteria are underlined. The gene(s) in parentheses are those that emerged from the more conservative analysis.

The general functional enrichment profiles of gF and gC were different. No function emerged as notably significant in gF, a finding which was further emphasized by the ‘flat’ multiple-testing corrected BH P-value of 0.18 for all significant annotations. The most significantly enriched functional annotation identified was ‘synaptic fatigue of synapse’ (FET P = 0.0021), followed by ‘apoptosis of spinal cord cells’ (FET P = 0.0041) and ‘inhibition of neurons’ (FET P = 0.0053) (Table3). However, the predominance of function annotations relating to ‘quantity’, particularly of neurons and of other structures of the nervous system, was notable (Table3; Fig. S2; Dataset S3). Other functions were related to the quantity (e.g. formation, loss, survival and apoptosis), quality (e.g. synaptic fatigue, degradation, atrophy and myelination) or morphology of neurons or related structures. Long-term potentiation (LTP) was another repeated function for gF (Fig. S2). One gene, brain-derived neurotrophic factor (BDNF), which showed gene-based association with gF (Dataset S2), contributed to the enrichment of several significant function annotations for gF (Table3 and Dataset S3). As BDNF is one of the most extensively studied genes in the field (Green ), it could bias the IPA results, which are based on peer-reviewed publications. However, when IPA was run without BDNF, most of the significant findings withstood its exclusion (Table3 and Dataset S3), and the general functional enrichment profile observed in gF was preserved. A manual examination of the chromosomal positions of the sets of genes leading to the enrichment of the significant annotations listed in Tables3 ensured that these significant enrichments are due to independent association signals and not due to the same association signal represented by clusters of genes in LD.

Table 3

Top 25 IPA FAs for gF

Function	FA	FET P-value	BH P-value	No. of genes	Genes
Synaptic fatigue	Synaptic fatigue of synapse	0.0021	0.18	2	BDNF, SYN3
Apoptosis	Apoptosis of spinal cord cells	0.0041	0.18	4	BAX, BCL2, BDNF, FAM134B
Inhibition	Inhibition of neurons	0.0053	0.18	4	BDNF, NPPA, NPPB, NRG1
Development	Development of sensory nervous system	0.0062	0.18	2	BDNF, KIF1A
Inhibition	Inhibition of pyramidal neurons	0.0062	0.18	2	BDNF, NRG1
Concentration	Concentration of arachidonic acid	0.0062	0.18	2	DGKE, KNG1
Endocytosis	Endocytosis of synaptic vesicles	0.0067	0.18	4	AMPH, CABIN1, ITSN1, SNCA
Damage	Damage of cortical neurons	0.0067	0.18	3	BAX, BDNF, GRIN3A
Quantity	Quantity of central nervous system cells	0.0090	0.18	9	ACHE, AVPR1B, CCND2, EGFR, GDA, IGFBP6, LEF1, LHX5, TSC1
LTP	LTP	0.0092	0.18	13	B3GAT1, BDNF, CDH1, CHRNA7, CYP19A1, DGKE, LRP1, LRP8, NRG1, PLG, SNCA, VAV2, VAV3
Apoptosis	Apoptosis of dorsal root ganglion cells	0.0093	0.18	3	BAX, BCL2, FAM134B
Mobilization	Mobilization of Ca2+	0.0093	0.18	3	BDNF, PROK2, TRPV1
Morphogenesis	Morphogenesis of neurites	0.0097	0.18	22	ACHE, ATXN10, BDNF, BSN, CNTN4, CNTNAP2, CTNNA2, CTNND2, CYP19A1, EGFR, GDA, IGF1R, MAPK8, NRG1, PRKG1, PTPRM, RND1, SGK1, SLC18A3, SLIT1, SYNE1, UHMK1, (TLR7)
Growth	Growth of dendrites	0.0099	0.18	5	BDNF, CTNND2, CYP19A1, NRG1, SLIT1
Morphology	Morphology of dendrites	0.0099	0.18	5	BDNF, DISC1, GRIN2D, NPAS3, NRG1
Quantity	Quantity of acetylcholine	0.012	0.18	2	ACHE, SLC18A3
Survival	Survival of dorsal root ganglion	0.012	0.18	2	BAX, BDNF
Synaptogenesis	Synaptogenesis of brain cells	0.012	0.18	2	BDNF, CYP19A1
LTP	LTP of granule cells	0.012	0.18	2	BDNF, DGKE
Sensitization	Sensitization of neurons	0.012	0.18	2	BDNF, KNG1
Migration	Migration of GABAergic neurons	0.012	0.18	2	GFRA1, SLIT1
Quantity	Quantity of astrocytes	0.012	0.18	4	ACHE, EGFR, IGFBP6, TSC1
Cell viability	Cell viability of motor neurons	0.015	0.18	5	BAX, BCL2, BDNF, GFRA1, REG3G
Loss	Loss of motor neurons	0.016	0.18	3	BCL2, BDNF, GFRA1
Quantity	Quantity of nerve ending	0.016	0.18	3	BDNF, SLC18A3, SNCA

For each general function, the specific FA is shown together with the corresponding FET P-value, the BH multiple-testing corrected P-value, the number of genes and the symbols of the genes responsible for the enrichment signal in the FA. FAs that remained significant after excluding BDNF are in bold. The functions and FAs that remained significant in the IPA analysis of the genes that passed the more stringent criteria are underlined. The gene(s) in parentheses are those that emerged from the more conservative analysis.

In order to gauge the robustness of the above findings in the absence of the gene-length bias that may have been introduced with the relaxed replication criteria, the CNS-specific IPA analysis was also performed on gene lists constructed on the basis of more stringent replication criteria. For gF, the single-marker analysis identified two genes which contained replicated SNPs that met the meta-analysis P ≤ 1 × 10−5 threshold (Dataset S1). The gene-based analysis identified 178 genes with gene-based P ≤ 0.05 using either the modified Sidak approach or Brown's approximation method (Dataset S2). Together, this resulted in a total of 180 genes available for IPA analysis for gF. For gC, 5 and 224 genes were identified via the single-marker and gene-based analyses, respectively, resulting in 225 unique genes available for IPA analysis. Nine genes were common to both gF and gC.

The IPA results emerging from the abridged gene lists were not as impressive in terms of the number of significant unique function annotations (N = 26 and 8 for gF and gC, respectively) and in terms of their general enrichment profiles for gF and gC (Dataset 3). Also, for both traits, different function annotations emerged as most significant, namely ‘neuritogenesis’ for gF (FET P = 2.8 × 10−4) and ‘hypoplasia of cerebellar vermis’ (FET P = 0.0015) (Dataset 3). However, the main functions that were found to distinguish gF from gC in the original, less conservative analysis remained significant. ‘LTD’ and ‘synaptic depression’ remained unique to gC albeit at a reduced significance level (FET P = 0.018 and 0.047, respectively) (Table2 and Dataset 3). For gF, functions which related to quantity, quality or morphology of neurons and synapses still predominated. These included ‘quantity’, ‘morphogenesis’ and ‘development’ of neurons, synapses and other brain structures (Dataset 3). No function annotations were common to both traits in this analysis.

Finally, returning to the original gene lists based on the relaxed criteria, IPA was run again on both traits, using information from all cell lines, with the exception of cancer cell lines, in the Ingenuity® Knowledge Base. This helped to evaluate the primary CNS-restricted results relative to more global annotation, to identify annotations that may be related to cognitive ability and to enable the inclusion of published studies of brain-related traits that were executed in non-CNS cell lines, such as lymphoblastoid cell lines (Gladkevich ). The significance of the original, CNS-restricted findings was largely unaffected. The most significant annotations for gF in this analysis fell under the disease and disorder category of ‘Cardiovascular Disease’ (e.g. ‘vascular disease’, FET P = 7.1 × 10−6; Dataset S3). In addition, the function annotations ‘development of brain’, ‘behavior’, ‘cognition disorders’, ‘cognitive impairment’ and ‘schizophrenia’ were also significantly enriched in gF (Dataset S3). For gC, ‘synaptic depression’ remained the most significant annotation. However, it was followed by the ‘Cardiovascular Disease’ function annotation ‘coronary artery disease’ (FET P = 1.8 × 10−6; Dataset S3). The significance of the enrichment of ‘schizophrenia’ also marginally improved (from FET P = 0.042 to FET P = 0.026, Dataset S3).

Discussion

Novel polygenic approaches to analyzing GWAS data have greatly improved our understanding of complex traits and have captured more of the ‘hidden heritability’ (Davies ; Maher 2008; Manolio ; Visscher ) but have failed to identify the markers or genes implicated. Pathway-based techniques are complementary polygenic methods that support biological analysis of GWAS data. Gene-based methods, which are ideal for pathway-based approaches (Liu ), permit the aggregation of SNPs of smaller effect and test the gene as a whole, potentially increasing power (Neale & Sham 2004). These analyses are also more permissive to locus heterogeneity, where multiple independent variants within a locus have independent effects on a trait (Christoforou ; Yang ), and they use prior biological knowledge, facilitating a more meaningful interpretation of data (Wang ). Therefore, if the genetic signals underlying gF and gC, two highly heritable and polygenic traits, cluster in known biological pathways, it should be possible to use pathway-based analyses to identify which biological processes are most strongly implicated in these subcomponents of g. We particularly wanted to determine whether gF and gC could be differentiated at the level of biological pathways or functions as they have been shown to be distinct in terms of development and brain structure. Taking a ‘bottom-up’ approach (Liu ), we analyzed two GWAS datasets using both single-marker and gene-based analyses to construct gene lists for IPA analysis. Although some significantly enriched functions overlapped in the two traits, the overall picture suggested distinct functional enrichment profiles, supporting the gF–gC distinction at the genetic level.

The most statistically significant finding was the enrichment in gC of genes involved in synaptic depression and LTD. This finding, which was unique to gC, also survived the FDR multiple-testing correction at the more conservative 0.05 threshold and the application of more stringent replication criteria for the construction of gene lists. It was also the top finding when the pathway analysis included non-CNS-related tissue types and conditions. Long-term depression and LTP are the major forms of long-lasting synaptic change in the mammalian brain (Collingridge ). Long-term depression is involved in synaptic pruning during development (Peineau ) and is thus important in adult neuroplasticity. Synaptic LTD is mediated by the effect of l-glutamate and other neurotransmitters on several types of receptors. The resulting synaptic plasticity is necessary for hippocampus-dependent learning and memory, certain types of behavioral flexibility and novelty detection (Collingridge ). Long-term depression deficits have been associated with reduced working memory and reversal of memory performance in rats. Studies focusing on medial temporal lobe regions have revealed a role of LTD in memory stabilization (in the amygdala) (Migues ) and recognition memory (in the perirhinal cortex) (Winters & Bussey 2005), implicating LTD in the development and maintenance of knowledge representations.

For gF, the association signal was less striking in terms of statistical significance. gF was predominantly characterized by genes that control the quantity, morphology and integrity of neurons and synapses. These factors affect the quality and efficiency of neuronal signaling (Brown ), which, in turn, affect cognition, as indicated by reduced activation in individuals with higher cognitive ability (Prat ; Reichle ). As gF declines in cognitive aging, the enrichment of these functions is consistent with the decline in synapse number, brain volume and white matter integrity in the aging brain (Fjell ; Morrison & Hof 1997; Walhovd ; Westlye ).

Only four function annotations were common to both traits. These were related to development and structural aspects (e.g. microtubule dynamics and dendrite formation), which play central roles in synapse formation and are thus likely to be important for cognitive function (Bramham ). The genes leading to enrichment of three of these functions were different for gF and gC, reflecting the gain of information that gene set or pathway-based approaches offer by enabling the identification of overlap between related traits or replication of the same trait.

The function annotation ‘schizophrenia’ was also significantly enriched in gF and gC, indicating that the identified genes have also been implicated in schizophrenia by other studies. A polygenic risk score for schizophrenia was recently shown to be associated with lower IQ at age 70 and greater decline in IQ level in one of the CAGES cohorts (Mcintosh ), confirming previous observations at the single candidate level of a genetic relation between general cognition and schizophrenia (Toulopoulou ). Furthermore, enrichment of the terms ‘behavior’, ‘cognition disorders’ and ‘cognitive impairment’ suggests that the genes identified here have been implicated in other studies of cognition. Finally, annotations relating to ‘Cardiovascular Disease’ ranked in the top two in both gF and gC. Cognitive dysfunction is well documented in patients with cardiovascular disease (Vogels ).

The heritability for gF and gC are similar (Davies ), but the strength of association, extent of replication and strength of enrichment were all greater for gC than for gF. There are several possible reasons for this. First, gF may be more heterogeneous at the phenotype level. While gC was defined as the standardized score of a single measure in each sample, gF estimation was based on a hierarchical principal component analysis from an array of specific measures. The idea that the psychometric structure of gF is more heterogeneous is consistent with another major model of intelligence (Johnson & Bouchard 2005; Vernon 1964), in which verbal abilities are retained as a single second stratum component, whereas fluid abilities are separated into perceptual and mental rotation skills, thus reflecting higher phenotypic complexity. Second, assessment of gF was not based on identical subsets of specific tests in NCNG and each of the CAGES subsamples. However, the different batteries of cognitive tests yield almost identical estimates of general intelligence (Johnson ). Third, the age differences between the samples could have a differential effect, given that gF decreases with age while gC is relatively stable, even though the correlation between them increases with age (Li ). Finally, the difference in enrichment profiles may be magnified by gene-length bias, given that a stronger association was identified for gF when IPA was re-run using the genes that passed the more stringent criteria.

The main complicating factor in this and other studies of brain-related traits is gene-length bias, as it presents an analytical ‘Catch-22’. While it is important to correct genes for their length, or equivalently, the number of SNPs tested, to control for false findings, doing so risks over-penalizing and thus eliminating the most relevant candidate genes and therefore pathways for intelligence. It is recognized that gene length is related to functional complexity (Xu ), and it is known that brain-expressed genes involved in relevant neuronal processes and/or associated with autism and intellectual disability are substantially longer (King ; Michaelson ). Also, it has been suggested that longer genes are larger targets and therefore more prone to random mutation and are subject to different transcriptional mechanisms that may be functionally relevant to brain-related traits (King ; Solier ). Thus, in order to balance the competing Type I and II error rates, we focused our analysis on a well-accepted and recommended discovery-replication approach (Jia ), running the replication in one of the largest samples of its kind (the CAGES). Strict replication criteria were used in the single-marker analyses (Konig 2011). For the gene-based analyses, a robust permutation-based approach, which accounts for LD structure and SNP number, was used to select genes for replication, the criteria for which were subsequently relaxed to avoid over-penalizing the larger, more relevant genes. In order to then assess the robustness of these findings, IPA was also run on genes that survived the more stringent approach that corrects for SNP number. As the power of the IPA analysis was greatly compromised by the ∼sixfold reduction in the number of genes available for the analysis, it was not surprising that the evidence for the genetic distinctiveness of the two traits was weaker. Nevertheless, the main distinguishing features between gF and gC were upheld. This suggests that the true genetic architecture of these two traits may lie in the middle, and is inaccessible using current approaches which either ignore or over-penalize for the length of the gene.

As larger consortia form and pathway analyses continue to improve, a better understanding of the genetic architecture of gF and gC will emerge (Khatri ; Lencz ). Meanwhile, our study serves as a starting point, supporting the gF–gC distinction at the genetic level and critically converging with the findings of developmental and cognitive neuroscience studies. The specific function annotations, or pathways, identified are worthy of further replication and interrogation, using, for example, ‘top-down’ approaches such as gene set enrichment analysis (Fernandes ) to test specific gene sets constructed on the basis of these findings.