Hierarchical Classes Analysis (HICLAS): A novel data reduction method to examine associations between biallelic SNPs and perceptual organization phenotypes in schizophrenia.

The power of SNP association studies to detect valid relationships with clinical phenotypes in schizophrenia is largely limited by the number of SNPs selected and non-specificity of phenotypes. To address this, we first assessed performance on two visual perceptual organization tasks designed to avoid many generalized deficit confounds, Kanizsa shape perception and contour integration, in a schizophrenia patient sample. Then, to reduce the total number of candidate SNPs analyzed in association with perceptual organization phenotypes, we employed a two-stage strategy: first a priori SNPs from three candidate genes were selected (GAD1, NRG1 and DTNBP1); then a Hierarchical Classes Analysis (HICLAS) was performed to reduce the total number of SNPs, based on statistically related SNP clusters. HICLAS reduced the total number of candidate SNPs for subsequent phenotype association analyses from 6 to 3. MANCOVAs indicated that rs10503929 and rs1978340 were associated with the Kanizsa shape perception filling in metric but not the global shape detection metric. rs10503929 was also associated with altered contour integration performance. SNPs not selected by the HICLAS model were unrelated to perceptual phenotype indices. While the contribution of candidate SNPs to perceptual impairments requires further clarification, this study reports the first application of HICLAS as a hypothesis-independent mathematical method for SNP data reduction. HICLAS may be useful for future larger scale genotype-phenotype association studies.

Introduction

A significant limitation of single nucleotide polymorphism (SNP) association studies is that testing efficiency is affected by the number of SNPs analyzed and sample size. SNP association is also influenced by genotype and risk allele frequency (Bhangale et al., 2008, Zondervan and Cardon, 2004). Selection of the most informative SNPs may help maximize power of common variants in association with phenotypes (Hinds et al., 2005). However, schizophrenia is characterized by significant genetic heterogeneity (Hallmayer et al., 2005, Owen et al., 2005, Sebat et al., 2009). Therefore, selection of individual SNPs based on hypotheses alone may not capture considerable genotypic variation associated with different study populations.

Hierarchical Classes Analysis (HICLAS) is a method for representing set-theoretical patterns for two-way, two-mode binary matrices (De Boeck and Rosenberg, 1988). HICLAS was employed in this study to represent SNP (allele) covariation patterns, and on the basis of those patterns, reduced the number of SNPs analyzed in relation to specific phenotypes. HICLAS assumes that SNP allelic distributions are modeled in a binary array. That is, for any given biallelic SNP (e.g. rs3924999), patient genotypes were represented as 00, 10 (01), or 11. 00 denotes an individual’s genotype is homozygous for the first allele (e.g. TT), 10 (or 01) denotes the individual’s genotype is heterozygous (e.g. TC), and 11 denotes the genotype is homozygous for the second allele (e.g. CC). To model a patient with respect to 6 SNP genotypes (12 alleles) there would be 12 binary entries, such as ‘001111001010’. The representation of the genotype (allele) patterns in a sample of 90 patients is accomplished by concatenating the data from individual patients row-wise, yielding a 90 × 12 binary matrix.

However, with rare exceptions, multiple algebraic decompositions exist for a two-way binary matrix. Unlike Boolean factor analysis (also compatible with two-way binary matrices), the HICLAS model has compatible and incompatible decompositions. The devised algorithm guarantees that the decomposition is compatible to the set-theoretical formulation of the model (De Boeck and Rosenberg, 1988). The equivalence and order relations (superset-subset) among object classes (here patients) and attribute classes (here SNPs) are summarized, as well as their association to each other. A HICLAS model of a two-way binary array is compared to the actual data array using a Jaccard measure (goodness of fit of the model to the data), ranging from 0 to 1, where 1 indicates perfect fit. Previous simulations have shown that HICLAS is able to recover the a priori categorical structure of the data matrix when random error is added (De Boeck and Rosenberg, 1988). Thus, HICLAS can identify patterns among row (patients) and column (allele) bundles, making it an optimal SNP/allele data reduction strategy for subsequent analyses in relation to various phenotypes.

In addition to significant genetic heterogeneity (Hallmayer et al., 2005, Sebat et al., 2009), schizophrenia has considerable phenotypic heterogeneity as specific perceptual, cognitive and behavioral abnormalities are only observed in patient subpopulations (Carpenter and Buchanan, 1994, Heinrichs, 2001, Raffard and Bayard, 2012). While the assessment of many perceptual and cognitive domains in schizophrenia is susceptible to generalized deficit confounds (Carter, 2005, Chapman and Chapman, 1978, Knight and Silverstein, 2001, MacDonald and Carter, 2002, Silverstein, 2008), several specific visual processes have been assessed in a manner that avoids many confounds (Dakin et al., 2005, Dima et al., 2009, Keane et al., 2013a, Keane et al., 2013b). One example is visual perceptual organization (Silverstein and Keane, 2011, Silverstein et al., 2013, Uhlhaas and Silverstein, 2005).

Perceptual organization refers to the binding of individual stimulus features into lines, edges, surfaces, and object representations (Place and Gilmore, 1980, Silverstein et al., 1996, Silverstein et al., 2000). Importantly, perceptual organization impairments are observed in schizophrenia independent of medication effects (Silverstein and Keane, 2011), and can be revealed as superior performance to control groups in specific psychophysical paradigms where prepotent grouping of targets and distractors interferes with the performance of healthy subjects (Knight and Silverstein, 2001, Place and Gilmore, 1980). Perceptual organization impairments in schizophrenia are consistently observed and associated with poor premorbid functioning, treatment response, and functional outcomes (Silverstein et al., 1998, Silverstein et al., 2000, Uhlhaas and Silverstein, 2005), suggesting that they may represent a severe illness subtype biomarker for schizophrenia (Farmer et al., 1983, Sham et al., 1996, Wickham et al., 2001).

Because perceptual organization deficits in schizophrenia are associated with impairments in cognitive organization (i.e., thought disorder, inappropriate affect, etc.), it is hypothesized that they reflect an aspect of a widespread reduction in cognitive coordination – or the ability to modulate signal processing based on current spatial and/or temporal contexts – in schizophrenia (Phillips and Silverstein, 2003, Phillips et al., 2015). Animal and healthy human studies of perceptual organization indicate that it is subserved by neural synchrony (Uhlhaas, 2013, Uhlhaas and Singer, 2006). Neural synchrony modulates spatial and temporal integration in cognitive processing, (Uhlhaas and Singer, 2010) and relies on NMDA and GABAergic functioning (Bartos et al., 2007, Phillips et al., 2015, Phillips and Silverstein, 2003, Phillips and Silverstein, 2013, Uhlhaas and Silverstein, 2005, Silverstein and Keane, 2011). Moreover, NMDA and GABAergic circuits are dysregulated in schizophrenia (Lewis and Moghaddam, 2006, Lisman et al., 2008, Moghaddam, 2003, Poels et al., 2014). Specifically, loss of parvalbumin positive GABAergic interneurons has been found to reduce neural oscillations (Lodge et al., 2009, Spencer, 2009, Woo et al., 2010) leading to cognitive symptoms, (Cho et al., 2006) including perceptual deficits (Uhlhaas and Singer, 2010, Uhlhaas et al., 2006a). Genetic variations in these pathways have been also observed in schizophrenia association studies (Cherlyn et al., 2010, Petryshen et al., 2005).

Therefore, three candidate genes were given precedence in relation to the hypothesized perceptual organization neural circuitry. Glutamate decarboxylase 1 (GAD1), a regulator of GAD67 (GABA synthesizing enzyme), was selected. GAD67 protein expression has been shown to be reduced in postmortem schizophrenia brains (Addington et al., 2005, Guidotti et al., 2000) especially in cortical areas (Coyle, 2006, Lewis et al., 1999). Lower GABA concentrations have also been observed in the visual cortex in schizophrenia (Yoon et al., 2010). Neuregulin 1 (NRG1), was selected since it has been associated with schizophrenia in multiple studies (Stefansson et al., 2002, Stefansson et al., 2003). Nrg1/ErbB4 signaling regulates GABAergic transmission in the adult cerebral cortex, subsequently influencing inhibitory cortical function (Rico and Marín, 2011) and regulation of NMDA receptors in the prefrontal cortex (Gu et al., 2005). Dystrobrevin binding protein 1 (DTNBP1), was selected since dystrobrevin regulates expression of the NMDA NR2A receptor subunit in hippocampal and cortical regions (Blake et al., 1998, Tang et al., 2009).

The goal of this study was to establish if HICLAS could be employed as a novel SNP data reduction application for determining the strength of links between SNPs in glutamatergic and GABAergic pathway genes and perceptual organization deficits in schizophrenia. A priori SNPs from 3 candidate genes (NRG1, DTNBP1 and GAD1) were analyzed using HICLAS to determine if HICLAS selected SNPs were associated with incidence of perceptual impairments.

Materials and methods

Study participants

The study was approved by the Rutgers–Robert Wood Johnson Medical School Institutional Review Board. All study participants provided written informed consent. The recruitment, diagnostic and inclusion/exclusion procedures were previously reported (Joseph et al., 2013). Briefly, the Diagnostic Interview for Genetic Studies (DIGS) (Nurnberger et al., 1994) was administered and medical records were reviewed to determine if participants met DSM IV-TR (APA, 2000) criteria for schizophrenia or schizoaffective disorder. Participants with current substance use, mental retardation, neurological disorders, other primary psychiatric diagnoses or poor performance on attentional-control stimuli (see Section 2.7 Perceptual task data analyses) were excluded. The demographic and clinical composition of the HICLAS selected sample (40 African American and 50 Caucasian patients) is shown in Table 1.

Table 1

Demographic and clinical composition of the study sample.

Factor (N = 90)	M	SD
Age (Years)	47.0	10.9
Age Psychosis Onset (Years)	22.4	7.6
Total Chlorpromazine Equivalent (mg)	542.1	474.9
Participant Education Level (Years)	12.5	2.5
Mother Education Level (Years)	12.4	2.7
Father Education Level (Years)	13.4	3.1
Shipley Vocabulary Subtest Score	90.2	13.4
PANSS Positive Factor Score	11.4	2.8
PANSS Negative Factor Score	16.3	4.3
PANSS Depression Factor Score	14.3	4.0
PANSS Cognitive Factor Score	13.9	3.4
PANSS Disorganized Factor Score	7.5	2.3
PANSS Excitement Factor Score	10.0	2.2
PAS Overall Score	2.7	.81
PAS Social Sexual Factor Score	4.9	3.7
Estimated Visual Acuity Both Eyes	20/32	–
Sex (% Female)	37.8	–
Race (% African American)	44.4	–
Handedness (% Left Handed)	14.4	–
Schizoaffective (%)	37.8	–
Smoking Status (% Current Smoker)	48.9	–
Visual Hallucinations (% Current)	14.1	–
Outpatient/Partial/Acute Program (%)	50.0/31.1/18.9	–

Note: PANSS = Positive and Negative Syndrome Scale, PAS = Premorbid Adjustment Scale.

Clinical assessments

Symptoms occurring 2 weeks prior to study enrollment were assessed with the Positive and Negative Syndrome Scale (PANSS) (Kay et al., 1986). Psychosocial development was evaluated using the Premorbid Adjustment Scale (PAS) (Cannon-Spoor et al., 1982).

Study participants had stable antipsychotic medication dosages which were converted to chlorpromazine equivalents based on published standards (Andreasen et al., 2010). The Shipley Institute of Living Scale vocabulary subtest (Zachary, 1991) was used to generate full-scale premorbid IQ estimates. Current visual acuity was estimated with a Snellen chart.

DNA extraction and genotyping

Eight SNPs from 3 candidate genes (GAD1, NRG1 and DTNBP1) were selected based on the following criteria: 1) SNPs had previously shown biallelic variation; 2) SNPs were previously reported to have a Minor Allele Frequency (MAF) of 5% or greater in Caucasians and were targeted by the HapMap project; 3) SNPs had previously shown positive association to schizophrenia spectrum disorders; 4) SNPs were within genes previously linked to NMDA or GABAergic neurotransmitter systems; 5) SNPs were also considered based on the SZgene www.szgene.org/ meta-analysis of genetic studies for schizophrenia (Allen et al., 2008).

Candidate SNPs were then uploaded to Illumina's Assay Design Tool (Illumina, San Diego, California, USA) http://www.illumina.com/ for probe and panel design (GS0013878-OPA). SNP compatibility was based on design score, design rank, MAF and validation status. SNPs were verified using the dbSNP [Bethesda (MD): National Center for Biotechnology Information, National Library of Medicine dbSNP database (version 131, February 2010)] http://www.ncbi.nlm.nih.gov/SNP/.

Saliva samples were collected from study participants using Oragene kits (DNA Genotek Inc., Ontario, Canada) and sent to the Toronto Center for Applied Genomics (TCAG) for DNA extraction and genotyping. TCAG was blind to participant diagnosis and all phenotypic data. An Autopure LS Gentra/Qiagen DNA extractor running Puregene chemistry was used to extract the DNA which was hydrated in 10 mM Tris-HCL pH 8.0, 1 mM EDTA. DNA concentration was quantified using a flurometer and Hoescht dye.

The extracted DNA samples were processed in 96-well plates with 4 genotyping control samples per plate. 250 ng of genomic DNA underwent allele specific oligonucleotide hybridization followed by extension and ligation. A universal PCR (primers labeled with Cy2 or Cy3) step for the 8 loci followed. The amplified products were hybridized onto a GoldenGate® Genotyping Universal BeadChip and scanned using Illumina iScan according to the manufacturer’s protocol (Fan et al., 2006).

Perceptual tasks

The perceptual stimuli employed, and corresponding psychophysical derivation of performance metrics, have previously been described, and are summarized below (Joseph et al., 2013, Keane et al., 2014). Both tasks have demonstrated good internal consistency, test–retest reliability and minimal practice effects (Pennefather et al., 1999, Silverstein et al., 2012, Strauss et al., 2013). Experimental stimuli were presented on LED monitors (60 Hz) at three testing sites. The viewing distance ranged from 620 to 650 mm at each site so that individual pixels subtended .025° of visual angle square. Stimuli were displayed at (achromatic) intensities of 59 cd/m2 (black) or 76 cd/m2 (white), as verified with a Konica Minolta LS-100 luminance meter.

Contour integration

The Jittered Orientation Visual Integration (JOVI) task is a test of contour integration that determines ability to integrate Gabor elements into a perceptual whole (Silverstein et al., 2012). Participants were shown static Gabor elements forming an oblong shaped contour embedded in a display of randomly oriented Gabor elements. Perceptual organization was manipulated by adding orientation jitter to the Gabor elements forming the contours, across 6 levels: +/− 0°, 7–8°, 9–10°, 11–12°, 13–14°, and 15–16°. For all stimuli, the ratio of the density of adjacent background elements to the density of adjacent contour elements was 0.9. At this level, adjacent contour elements are farther apart than adjacent background elements, and thus contour identification cannot be accomplished via detection of density cues, and requires perceptual organization.

The JOVI is a symmetric 1 alternative forced choice task in which participants responded whether the narrow end of the oblong contour was pointing left or right for each trial (Fig. 1). Each stimulus was presented for 2 s followed by a 1 s inter stimulus interval during which responses were no longer recorded. 48 stimulus trials per jitter condition were presented in blocks of 12 trials. Two types of catch stimuli (i.e., no errors expected) using 0° jitter were administered during each block to assess momentary attention lapses. One had curved lines drawn through the contours to highlight contour salience, and the other contained contour elements without background elements to eliminate distractor noise effects. The task stimuli were created using E-prime (Psychology Software Tools, Pittsburgh, PA).

Fig. 1

JOVI task stimuli. The top left panel of the figure is an example of a lower jitter degree condition presented to participants (7–8°). The top right panel of the figure shows the highest jitter degree presented (15–16°). The bottom left and right panels represent the catch trial stimuli included in each trial block to account for momentary attentional lapses.

Kanizsa shape perception

Stimuli consisted of four white sectored circles (diameter = 3.0°; wedge = 45°) centered at the vertices of an invisible square (side = 9.0°), which itself was centered on the screen (Fig. 2). The unrotated pac-men in the illusory condition formed a square, one third of which was physically specified (support ratio = .33) (Kellman and Shipley, 1991). Certain trials contained distractor lines (dimensions = 4.0 × 0.1°), which were centered between the sectored circles and had a length equal to 2/3 of the illusory edge. A fixation point appeared at the screen center at the beginning of each trial.

Fig. 2

Kanizsa shape perception stimuli and trial sequence. (A) Participants discriminated illusory or fragmented squares, which were accompanied by distractor lines for some half of the trials. (B) The task was to say left/right for the fragmented condition or fat/thin for the illusory condition.

One half of the task consisted of the illusory condition, and the other half the fragmented condition (see Fig. 2A). The ordering of the conditions was counterbalanced across participants. In the illusory condition, the sectored circles were individually rotated clockwise or counter-clockwise by the same magnitude to form fat or thin shapes (Kellman and Shipley, 1991). For the fragmented trials, the elements were oriented downward (to prevent illusory contours) and were individually rotated to the right or left all in the same direction. A left/right task was chosen because it forced participants to make judgments on the lateral properties of the stimulus—similar to the illusory condition.

For each half of the task, there were 64 practice trials and 84 non-practice trials, the latter half of which presented distractor lines. This number of practice trials (Keane et al., 2012, Zhou et al., 2008) was selected to acclimate participants to the presentation times and orientation differences. The first non-practice trial for each condition was excluded for threshold estimation since these trials were often missed by observers. Participants received a brief break between blocks and preceding the distractor line trials.

The trial presentation sequence (Fig. 2B) was similar to earlier studies (Keane et al., 2012, Keane et al., 2014, Ringach and Shapley, 1996, Zhou et al., 2008) and consisted of a 1000 ms black screen, a 200 ms target presentation, a 50 ms uniform black screen, and 300 ms mask (to cap stimulus processing time). Another black screen would linger until a response, after which an auditory beep sounded for a correct answer. To reduce keyboard press errors, participants verbally responded “left”/“right” or “fat”/“thin” after each trial, with the experimenter subsequently entering the participant’s response.

Task instructions were presented before and after the practice trials. On one screen, luminance-defined lines were drawn on the borders of the illusory shape, so that participants clearly understood “fat” vs. “thin”. On subsequent screens, starkly different fat/thin shapes (rotation = 10°) were shown individually, side-by-side, and then in temporal succession (period = 2 s). During practice trials, the target presentation time and rotational magnitude decreased incrementally (3200 ms, 1600 ms, 800 ms, 400 ms, and 200 ms; 10, 8, 6, and 4°) to acclimate participants to subtle shape differences and brief stimulus presentation.

Task difficulty depended on rotational magnitude, with larger rotations making the alternatives easier to distinguish. A Bayesian adaptive “Psi” method (Kontsevich and Tyler, 1999) recommended a rotational magnitude for each trial, based on performance on previous trials, to minimize uncertainty of the slope and threshold estimates of the psychometric function. Rotational magnitude was expressed in log units given the decelerating function relating this quantity to proportion correct (Zhou et al., 2008). The algorithm assumed a log-Weibull (Gumbel) function (Prins and Kingdom, 2009).

where ψ is the proportion correct, x is the rotational magnitude, α is threshold, β is slope, γ is the guess rate (.5), and λ corresponds to the proportion of accidental responses (assumed to be .03) (Wichmann and Hill, 2001). Threshold establishes the position of the sigmoidal curve along the abscissa and corresponds to the rotational magnitude (in log degrees) needed for 79.7% accuracy. The Psi method was selected because it makes no assumption about slope – which can vary by condition – and because it provides an efficient means for estimating two parameter psychometric functions, (Klein, 2001) yielding a reliable threshold estimate (± 2 dB) with as few as 30 trials. Task stimuli were created in MATLAB with the Psychophysics Toolbox (Brainard, 1997).

Clinical data analyses

Spearman correlations were calculated to examine associations among premorbid functioning, symptoms, antipsychotic dosage, visual acuity and perceptual scores, and these relationships have previously been reported (Joseph et al., 2013, Keane et al., 2014). PANSS syndromes were analyzed based on a medication stable five factor model, (Lindenmayer et al., 1994) including positive, negative, cognitive, excitement, and depression factors. A separate disorganization factor (Cuesta and Peralta, 1995) was derived and included: poor attention, conceptual disorganization and inappropriate affect (the latter not being an original PANSS item). For the PAS, the social–sexual functioning factor and an overall mean score (Cannon-Spoor et al., 1982) were calculated.

SNP genotyping and HICLAS analysis

SNP cluster plots were manually inspected with GenomeStudio v.2011 to determine genotypes using default parameters. SNPs were called for GenCall scores > 0.25. rs11542313 had a poor call rate and rs16876589 lacked allelic variation in our study sample, reducing our total number of candidate SNPs from 8 to 6 (rs10503929, rs3924999, rs1978340, rs3213207, rs1047631 and rs1747054) for subsequent HICLAS. The sample included for HICLAS had a genotype rate of 98.7%. The subject by SNP data were then arranged in a 90 × 12 data array and were analyzed using the HICLAS algorithm separately for ranks 1 through 12 (rank refers to the Schein rank of any given HICLAS model of the binary array; a solution in rank 3 is roughly analogous to a 3 factor solution).

Perceptual task data analyses

For the contour integration task, the total score across all jitter conditions (excluding catch trials) was the performance index as this score has a higher test–retest reliability compared to psychometric function threshold values (Silverstein et al., 2012). Eleven participants were unable to identify the catch trials at a rate of at least 83.3% – a rate significantly better than chance based on total catch trial number – and were excluded prior to HICLAS.

For Kanizsa shape perception, two metrics were of interest. One was global shape integration, corresponding to how well participants distinguished Kanizsa shapes relative to featurally similar fragmented shapes (without distractor lines). A lower relative threshold in the illusory condition demonstrates an enhanced capacity to take advantage of the Gestalt layout of the stimulus. The second was how well participants fill-in illusory contours. Filling-in measured how much participants responded to seemingly irrelevant information (distractor lines) placed near the filled-in paths. Filling-in was operationalized on the basis of distractor effects: the more that distractor lines impaired discrimination in the illusory relative to the fragmented condition, the more that filling in was assumed to occur. This metric was chosen because others have shown that distractor lines near the edges of Kanizsa shapes worsen illusory shape perception, but have little effect when illusory contours are not perceived (Keane et al., 2012, Keane et al., 2014, Ringach and Shapley, 1996, Zhou et al., 2008). Therefore, based on prior studies from our lab and others, we propose that the Kanizsa shape perception task employed for this study is assessing two primary perceptual processes.

Results

Clinical and perceptual task correlations

Total JOVI scores were significantly correlated with increased conceptual disorganization and poor premorbid social sexual functioning, replicating previous findings (Joseph et al., 2013, Schenkel et al., 2005, Uhlhaas et al., 2006b). Age and sex were included as covariates for all ANCOVA/MANCOVA analyses since the contour interpolation global shape metric suggested a trend level correlation with participant age, and also because novel sex differences were observed in this sample (Joseph et al., 2013). No significant correlations between estimated visual acuity and perceptual indices were observed in this study.

HICLAS selected SNPs

The SNP functional classes, minor allele and genotype frequencies for Caucasian and African American participants are shown in Table 2. HICLAS selected SNPs were in Hardy Weinberg Equilibrium (p > .05) for both Caucasian and African American participant groups. SNP genotype frequencies varied based on race: rs10503929 (χ (2, N = 86) = 7.3, p = .026), rs3924999 (χ2 (2, N = 86) = 8.5, p = .014), rs1978340 (χ (2, N = 86) = 7.2, p = .028).

Table 2

Minor allele and genotype frequencies in study sample for HICLAS selected SNPs.

Chr	Gene/SNP	Functional Class	MAF African American	MAF Caucasian	Genotype Frequency African American	Genotype Frequency Caucasian
2q31	GAD1rs1978340	5’ Flanking	0.125	0.266	CC 0.750CT 0.200TT 0.050	CC 0.510CT 0.449TT 0.041
8p12	NRG1rs10503929	Missense Methionine (T) to Threonine (C)	0.013	0.132	TT 0.975TC 0.025	TT 0.766TC 0.224CC 0.020
8p12	NRG1rs3924999	Missense Arginine (C) to Glutamine (T)	0.138	0.357	CC 0.750TC 0.225TT 0.025	CC 0.470TC 0.346TT 0.184

HICLAS SNP data reduction

The HICLAS solution in rank 3 yielded a Jaccard Index (goodness of fit) of .866. This solution was chosen over rank 2 (Jaccard Index = .793) because fewer SNPs were discarded. Rank 3 was chosen over rank 4 (Jaccard Index = .915) because the participant clusters and fit in the rank 4 did not change appreciably from how they clustered in rank 3. In addition, increasing the rank leads to a monotonic increase in fit.

The rank 3 solution generated 7 clusters of participants based on the presence/absence of 3 alleles from SNPs rs3924999 (NRG1), rs10503929 (NRG1) and rs1978340 (GAD1), as shown in Table 3. The other SNPs, rs3213207 (DTNBP1), rs1047631 (DTNBP1) and rs1747054 (DTNBP1), were eliminated. Hence, HICLAS reduced the total number SNPs to be analyzed in relation to perceptual indices from 6 to 3.

Table 3

Bundle patterns for HICLAS selected SNPs.

HICLAS Cluster	#Patients in Cluster	rs1978340	rs10503929	rs3924999
001	3			allele #1
010	2	allele #1
011	3	allele #1		allele #1
100	36		allele #2
101	14		allele #2	allele #1
110	11	allele #1	allele #2
111	16	allele #1	allele #2	allele #1

Genotypes for these three SNPs were then analyzed in relation to the perceptual task indices in 3 sets of ANCOVAs, with each SNP genotype as the 3-level independent variable (see Table 2 for SNP genotypes). The covariates included age and sex, and the dependent (perceptual organization) variables were total JOVI score, and the Kanizsa shape perception filling in and global shape metrics as shown in Table 4. The rs10503929 CC genotype was associated with poorer contour integration performance and reduced Kanizsa shape perception filling in abilities (higher Kanizsa shape perception filling in and more negative difference score of threshold scores) compared to TC and TT genotypes. The rs19783420 TT genotype was also associated with poorer filling in abilities for Kanizsa shape perception, as compared to CT and CC genotypes.

Table 4

ANCOVAs of HICLAS selected SNP and perceptual task indices.

Independent Variables	Dependent Variables	Covariates	F(2,83)	pa	ηp2	Power
rs1978340	JOVI Task Total Score	Age	2.0	.140	.047	.405
	Contour Interpolation Global Shape	Sex	.635	.533	.016	.153
	Contour Interpolation Filling In		4.3	.017	.098	.733
rs1053929	JOVI Task Total Score	Age	4.7	.012	.103	.771
	Contour Interpolation Global Shape	Sex	.537	.586	.013	.136
	Contour Interpolation Filling In		3.7	.028	.086	.668
rs3924999	JOVI Task Total Score	Age	.579	.563	.014	.143
	Contour Interpolation Global Shape	Sex	1.6	.219	.038	.319
	Contour Interpolation Filling In		3.1	.052	.072	.577

Due to the exploratory nature of this study, all p levels reported are uncorrected.

The results of the MANCOVA combining the 3 HICLAS selected SNPs are shown in Table 5.

Table 5

MANCOVAs for HICLAS selected SNPs and perceptual task indices.

Independent Variables	Dependent Variables	Covariates	Wilks’ λ	F(4,156)	pa	ηp2	power
rs1978340	Total JOVI Task Score	Age	.848	3.4	.011	.079	.838
rs3924999	Contour Interpolation Filling In	Sex	.920	1.7	.163	.041	.500
rs1053929			.837	3.6	.007	.085	.869

Due to the exploratory nature of this study, all p levels reported are uncorrected.

In addition, a MANCOVA combining the 3 HICLAS eliminated SNPs is shown in Table 6. Only HICLAS selected SNPs had a significant association to perceptual organization phenotypes.

Table 6

MANCOVAs for HICLAS eliminated SNPs and perceptual task indices.

Independent Variables	Dependent Variables	Covariates	Wilks’ λ	F(4,156)	pa	ηp2	power
rs3213207	JOVI Task Total Score	Age	.993	.147	.964	.004	.080
rs1047631	Contour Interpolation Filling In	Sex	.975	.496	.739	.013	.166
rs17470454			.971	1.2	.319	.029	.248

Due to the exploratory nature of this study, all p levels reported are uncorrected.

Discussion

HICLAS was employed as a novel application of a structural binary array model to reduce the number of candidate SNPs analyzed in relation to perceptual organization phenotypes in a schizophrenia sample. HICLAS reduced the total SNPs analyzed in relation to phenotypes from 6 to 3. Although the sample was restricted to Caucasian and African American participants, a significant study limitation is the small sample size that minimized power and did not allow for separate analyses by race. Future studies examining HICLAS selected SNPs in larger samples are needed to validate genotype-perceptual phenotype relationships.

The perceptual tasks evaluated in our study were selected based on the assumption that they represent a stable patient phenotype and a recent paper supports this for contour integration task performance (Feigenson et al., 2014). This suggests that the study tasks may be suitable for SNP-phenotype association analyses. HICLAS selected rs10503929 and rs1978340 genotypes were associated with Kanizsa shape perception filling in but not with global shape performance. Keane et al. (2014) suggest that filling in is an earlier perceptual stage compared with global shape discrimination, linked to the higher-level impairment of conceptual disorganization, indicating a more severe illness course. rs10503929 genotypes were also associated with contour integration performance. Contour integration impairments in schizophrenia involve hypo- and hyper-activation in different frontal regions (Silverstein et al., 2009), suggesting that template matching and decision making process abnormalities also contribute to task performance, and this may explain the relationships between both tasks and rs10503929.

Since contour integration and Kanizsa shape perception have distinct developmental trajectories, (Csibra et al., 2000, Káldy and Kovács, 2003, Kovacs et al., 1999) future studies may consider inclusion of genes related to neurodevelopment such as Fragile X Mental Retardation Protein (FMRP). FMRP has a key role in neuroplasticity (Fernández et al., 2013) and is transcriptionally repressed in schizophrenia, neurodevelopmental and mood disorders (Abel and Zukin, 2008). Kelemen et al. (2013) reported that FMRP protein levels are associated with contrast sensitivity, perceptual integration, and motion perception. FMRP is also thought to interact with glutamatergic and GABAergic pathways (D’Hulst and Kooy, 2007). Inclusion of SNPs from these genes may enable comprehensive modeling of epistatic contributions to perceptual organization phenotypes.

Although HICLAS has not been applied previously to genetic data, it offers some advantages over other genotype–phenotype association methods. First, it is a mathematical model unlike most set based SNP association analyses, which are thought to primarily consider linkage disequilibrium between SNPs (Liu et al., 2010). In addition, gene pathway based associations are usually not amenable to the inclusion of covariates and are subject to permutation biases (Wang et al., 2010). Other binary array models such as Boolean factor analyses are limited in that they do not consider the set-theoretical structure of the data array (De Boeck and Rosenberg, 1988). Therefore, HICLAS may be a valuable variable reduction method for future large scale SNP phenotype association studies.

Role of Funding Source

This study was partially funded by an American Psychiatric Foundation award to Dr. Steven M. Silverstein.

Contributors

JJ assisted with experimental design, collected the data, conducted literature searches, undertook statistical analyses and wrote the first draft of the manuscript. MAG performed the HICLAS and contributed to data interpretation, experimental design and editing of manuscript drafts. SMS designed and conducted the study, and edited the manuscript drafts.

Conflicts of Interest

The authors declare no conflicts of interest.

Acknowledgements

We thank Dr. Thomas V. Papathomas, Dr. Angus W. MacDonald III, and Dr. Linda M. Brzustowicz for their helpful guidance during the dissertation that led directly to this paper. We also thank Dr. Brian P. Keane for programming the Kanizsa shape perception task. We thank Dr. Tara Paton, Dr. Christian Marshall, Dr. Pingzhou Hu, Carol Ann Ryan and the staff at TCAG for their genotyping services. We also thank Dr. Brinda K. Rana and Dr. Paul DeBoeck for helpful feedback on an earlier draft of the manuscript. We acknowledge the Rutgers UBHC staff and patients for their support.