Technical Note: Informatics

DNA Copy Number and Loss of Heterozygosity
Analysis Algorithms
Detection of copy-number variants and chromosomal aberrations in GenomeStudio® software.
Introduction
Illumina has developed several algorithms for detecting copy number
variants (CNVs) and other structural variants from microarray data
(Table 1). These algorithms are available as individual software plug-ins
for the GenomeStudio Genotyping Module and can be downloaded
from the Illumina support page: http://support.illumina.com/
downloads/genomestudio-2-0-plug-ins.html. The plug-ins are usedwithin the CNV Analysis workbench, and results can be visualized
within the GenomeStudio Full Data Table, in the Illumina Genome
Viewer (IGV), or in a CNV region display window. This technical note
describes the function of these algorithms and how they can be
employed to analyze a chromosomal region of interest.

CNV Region Report

Table 1: GenomeStudio Copy
Number Algorithms
Algorithm

Function

CNV Region
ReportGenerates three separate CNV reports

cnvPartition

Calculates copy numbers with confidence
scores and generates CNV regions

Homozygosity
Detector

Autobookmarks samples with extended
tracts of homozygosity (single-sample
analysis only)

LOH Score

Estimates the likelihood of a region
exhibiting LOH

CNV Region Report is a software plug-in for GenomeStudio thatgenerates three separate CNV reports.
• Standard Report—Lists each CNV and loss of heterozygosity(LOH) region for each selected sample.
• Allele-Specific Copy Number Report—Estimates the allelespecific copy number for each probe entry (e.g., A- or AAB). In the
output file, the CN_GTYPE column is calculated using the CNV
Value (as determined by CNVPartition), the B Allele Frequency
data, the GTYPE (genotype column), and the theoretical B Allele
Frequency normal distributions for each copy number.
• PLINK CNV Input Report—Creates input files for some of the CNV
features of the PLINK genome-wide association study (GWAS) and
CNV analysis application1.

cnvPartition
The goal of the cnvPartition algorithm is to identify regions of the
genome that are aberrant in copy number using two Infinium® assay
outputs: the log R ratio (LRR) and B allele frequency (BAF). Because
LRR is the logged ratio of observed probe intensity to expected
intensity, any deviations from zero in this metric are evidence for copy
number change. BAF is the proportion of hybridized sample thatcarries the B allele as designated by the Infinium assay. In a normal
sample, discrete BAFs of 0.0, 0.5, and 1.0 are expected for eachlocus (representing AA, AB, and BB).
Deviations from this expectation are indicative of aberrant copy
number. For example, if a locus has a BAF of 0.66, this might indicatethat there are two copies of the B allele and one copy of the A allele
2
present in the sample 2 + 1 ≈ 0.66 . Analyzing both of these metrics
provides stronger resolution for detecting true copy number changes.

Copy Number Estimation
cnvPartition models LRRs and BAFs for each of 14 different copy
number scenarios as simple bivariate Gaussian distributions (Table 2).
Modeling copy number in this way allows for computation of apreliminary copy number estimate for each assayed locus bycomparing its observed LRR and BAF to values predicted from eachof the fourteen models. Specifically, the likelihood of observing agiven LRR and BAF under each of the 14 models is calculated. For
example, to compute the likelihood of a particular LRR and BAF givena genotype of AAB (LAAB), the AAB parameters from Table 2 and thestandard normal density are used:
1( LRR − 0.3)2
1(BAF − 13 )
+
exp −
exp −
2(0.18 2 )0.18 2π0.03 2π2 (0.032 )2

LAAB =

Likelihoods are also computed for other model genotypes listed inTable 1 with the exception of the homozygous deletion (DD). For
homozygous deletions, a very low LRR is expected, but the BAF may
be any value between zero and one. Therefore, the likelihood of adouble deletion (LDD) is calculated by the equation:
LDD =1
2 × 2πexp −(LRR − ( −5))22 (22 )

Technical Note: Informatics
Table 2: Genotypes Modeledby cnvPartition
Si = X1 + … + X i ,1 ≤ i ≤ n
Genotype

CN

LRR
Mean

LRR
SD

BAF
Mean

BAF
SD

DD

0

-5

2

NA

NA

A

1

-0.45

0.18

0

0.03

B

1

-0.45

0.18

1

0.03

AA

2

0

0.18

0

0.03

AB

2

0

0.18

0.5

0.03

BB

2

0

0.18

1

0.03

AAA

3

0.3

0.18

0

0.03

AAB

3

0.3

0.18

1/3

0.03

ABB

3

0.3

0.18

2/3

0.03

BBB

3

0.3

0.18

1

0.03

AAAA

4

0.75

0.18

0

0.03

AAAB

4

0.75

0.18

0.25

0.03

ABBB

4

0.75

0.18

0.75

0.03

BBBB

4

0.75

0.18

1

0.03

Parameters for each of the fourteen genotypes considered by cnvPartition
are shown. BAFs are modeled as a uniform distribution between zero andone for homozygous deletions (DD). All other distributions are modeled with
Gaussian distributions with the given parameters. The genotype AABB is
not modeled since this would represent two independent duplication eventsand rarely occurs in nature. (CN = copy number, DD = double deletion, SD= standard deviation)

These likelihoods are then summarized by four composite copy
number likelihoods:
L0 = LDD
L1 = LA + LB
L2 = LAA + LAB + LBB
L3 = LAAA + LAAB + LABB + LBBB
L4 = LAAAA + LAAAB + LABBB + LBBBB

where Lk denotes the likelihood of copy number k for integer valuesof k and the likelihood of a genotype for non-numeric values of k. The
preliminary copy number estimate (X) is defined as the average of thefive modeled copy numbers, weighted by their respective likelihoods:
X=

L1 + 2L2 + 3L3 + 4 L4L0 + L1 + L2 + L3 + L4

Breakpoint Identification
Preliminary copy number estimates are the inputs to the core
partitioning algorithm. The goal of partitioning is to identify regions of
the genome where the values of X are consistently higher or lower than2, the expected value for a diploid sample. To find an aberrant region,
the algorithm orders the X values by their position along a chromosome
and searches for the indexes i and j such that the values Xi . . .Xj aremaximally different than those outside this region. Thus, the algorithm
seeks to maximize |Zij| over all i and j with i < j, defined by the equations:

Z ij =

1
1+( j − i) (n − j + i)

−1/2

×

( S j − Si ) ( Sn − S j + Si )−( j − i)(n − j + i)

where n is the number of loci assayed on the chromosome.
An exhaustive search through all pairs of i and j scales quadraticallywith n and is therefore an inefficient process for use with Illumina
whole-genome genotyping products2,3. To simplify the calculations
required, cnvPartition uses a sliding window strategy to maximize|Zij|, but where j = i + w, with w the defined window size. After the
optimal window size value is found, the algorithm attempts to extend
the window in both directions to maximize the value of |Zij| further.
As implemented, the algorithm repeats this procedure for w = 4, 8,
16, and 32 then reports the i and j corresponding to the maximal|Zij| found. When a maximally different segment is found, |Zij| iscompared to a pre-determined threshold (default is 6). If the thresholdis exceeded, the boundaries are noted and the algorithm is appliedrecursively to the regions between 1 and i, i+1 and j, and j+1 throughn. The threshold of 6 was chosen as a default because it minimizes
false positives, particularly for short aberrations.

Copy Number Assignment to Partitioned Regions
The partitioning procedure results in a set of putative breakpoints
scattered across the genome. The next step is to assign a copy
number for each region lying between two consecutive breakpoints.
To do this, L0, L1, L2, L3, and L4 for each locus within the region areused. For each putative copy number (0–4), the logarithms of all Lk foreach k are summed. The k with the highest sum is the copy number
assigned to this region. For regions with copy numbers other than 2,
the algorithm also assigns a confidence score for the copy numberthat is called. The confidence score is defined as the sum of all logged
likelihoods in the region for the assigned copy number minus the sumof all log L2 values for loci in the region.

Additional Usage Notes• Regions with copy number = 1 on the X or Y for males are filteredfrom the CNV results.
• Probes that are designated as Intensity Only are treated differently thannormal probes. The B Allele Frequency is ignored for these probes.
• Y probes are not considered for samples designated as female.

Homozygous Region Detection
cnvPartition also includes a homozygosity detection algorithm thatruns separately from the partitioning algorithm already described. Thisalgorithm only runs on copy number 2 regions by default. Therefore,
it is sometimes called a Copy Neutral LOH Detector. The logic issimilar to that used in the homozygosity detector autobookmarking
plug-in (see next section). However, additional logic has been added
to simplify usage for the end user. Instead of adjusting the ChiSquarethreshold as in the autobookmarking plug-in, the user can simply
adjust the MinHomozygousRegionSize configuration parameter. By
default, this is set to 10 Mb based on empirical testing.

Technical Note: Informatics
Experimental Features
Several experimental features are included in cnvPartition v3.2.0.
These features are all disabled by default but can be activated byadjusting the associated configuration files. Configuration details foreach feature are available in the cnvPartition user document5. Thesefeatures include:• LOH detection for the entire genome. Previous versions only
allowed LOH detection to be run on regions with a copy number of 2.
• Log R Ratio adjustment for the Y chromosome. Y SNPs are
clustered using only males, so the log R Ratio indicates that the
copy number = 2. The adjustment lowers all log R Ratios on Y.
If Y chromosome SNP clusters are already adjusted, additional
adjustment could provide inconsistent results.
• Support for highly amplified genomes. These are common withmany cancer samples. Log R Ratios are adjusted upward based on
calculation of average genomic ploidy.
• A GC wave adjustment based on linear regression of LRR vs. GC
content in probes. GC waves can be present when using incorrectly
quantified DNA in the Infinium Assay, or they might be present inregions of high or low GC content.

CNV Display in GenomeStudio
The copy number values are then used to create CNV regions andbookmarks in GenomeStudio for visualization of aberrant regions. The
cnvPartition algorithm incorporates two user-definable thresholds for
optimization of CNV detection. The Confidence Threshold allows users
to filter out CNV regions that have low confidence values. The default
of 35 was determined empirically using normal HapMap samples on
the Illumina Human1M BeadChip. The Probe Gap Size Thresholdallows users to filter out CNV regions that are in large probe gaps,
such as centromeres. The default of 1,000,000 (1 Mb) was determined
empirically to help prevent CNV regions from being falsely detectedacross centromeres and other large gaps.

Algorithms for Automated Bookmarking
GenomeStudio can use several automated bookmarking algorithms.
These plug-in algorithms automatically scan data for the presence ofstructural aberrations and CNVs. Each algorithm employs a different
strategy to search for aberrations and CNVs. They are intended toassist with visually categorizing various types of aberrations present
in samples of interest. Automated bookmarking algorithms can beused as data-mining tools for the discovery of new regions, or toverify known regions of interest. Bookmarks can be edited whether
they are created manually or generated using the autobookmarking
tool. All bookmarks can be exported to share with other users. The
Homozygosity Detector is an autobookmarking plug-in available for
GenomeStudio.

Homozygosity Detector
The homozygosity detector algorithm can be used to autobookmark
samples with extended tracts of homozygosity (single-sample
analysis only). Homozygosity tracts may result from inbreeding,large-scale gene conversion, uniparental disomy (UPD), or
chromosomal deletions. Other factors such as population history or
low recombination rates may also contribute to creating extended
regions of LOH. This algorithm uses SNP frequencies to calculate

the expectation that a single SNP is homozygous in a single sample.
The algorithm then calculates the X2 value of the observation ofzero heterozygotes in N SNPs versus the expected number ofheterozygotes given the frequencies of the SNPs. In this algorithm,
there is a fixed cutoff significance (X2 = 23.5), which corresponds to
50 contiguous SNPs, each with a minor allele frequency (MAF) of 0.2.
This algorithm requires that each LOH region contains at least
50 homozygous SNPs by default. All LOH regions with more than
50 SNPs and X2 > 23.5 are bookmarked.

Homozygosity Detector Algorithm Process
1. The X2 threshold is preset to 23.5, but users can change this value.
2. The minimum number of SNPs per region is preset to 50, butusers can change this value.
3. The allele frequencies are calculated for each SNP before analysis.
The allele frequencies are used to calculate the expectation thateach SNP is heterozygous in a single sample, assuming HardyWeinberg equilibrium.
4. Next, for each sample, all of the genotypes on each chromosome
are scanned, and all of the contiguous regions without
heterozygous genotypes are located.
5. For each of these regions, the expected number of heterozygotes
is calculated by the equation:N

E het = ∑ 2fi ( l − fi )i =1

where N is the number of homozygous SNPs and fi is the frequencyof either of the SNP alleles in the general population. The X2 value isgiven as:X2 =( Nhom E hom )2 ( Nhet − E het )2
+
E hom
E het

where Nhom and Nhet are the number of homozygous andheterozygous genotypes respectively, and Ehom and Ehet are theexpected number of homozygous and heterozygous genotypes.
By definition, there are no heterozygotes, so the X2 value can be
simplified to:X2 =NE hetN − E het

where N is the number of SNPs with genotype calls.
6. Each segment that is more significant than the predefined or usersupplied X2 threshold value and has more SNPs than the predefinedor user-supplied minimum number of SNPs is bookmarked.

LOH Score Statistical Algorithm
The LOH Score column plug-in reports the likelihood of loss ofheterozygosity (LOH) existing in a region of interest. The LOH Score
algorithm scans data sets to determine and identify the presence ofLOH. Variances in LOH score can be plotted in the Chromosome Heat
Map or in the Illumina Genome Viewer (IGV).

Technical Note: Informatics
Figure 1: LOH Score Equation

Perror = 0.001

N

Π
Π

i=1
LOHScore = Log
10N
i=1Pi(gt i|LOH) ={Pi(gt i|NoLOH) ={Pi(gt i|LOH)
Pi(gt i|NoLOH)

Perror, gti = AB
1 − Perror, gti = AAorBB

hetfreq, gti = AB
1 − hetfreq, gti = AAorBB

Perror = genotyping error rate
hetfreq = mean heterozygote frequency
gti = genotype of locus i

If a chromosomal region is lost and LOH is observed, the only
expected genotypes are AA and BB. In this case, AB would beobserved only as a result of genotyping error. If there is no LOH, allthree genotypes are possible.
The LOH score is a measure of the likelihood that a SNP is exhibiting
LOH around a window over all N SNPs, where N is the number of
SNPs in a user-designated window size centered at the chromosomal
position of the SNP. The equation used to determine the LOH Score isshown in Figure 1. The recommended window size depends upon thedensity of probes on the product in use, as defined by the following
equation
Factor Number
Number of Probes
where Factor Number = 88,000

Minimum Window Size =The window size for a specific workspace may require optimization
depending upon the type of aberration under examination and thequality of the data.
If the number of heterozygote SNPs matches the prediction in thespecified window, the LOH score is 0. The LOH score increases if
there is an unexpectedly low number of heterozygote SNPs in thewindow. Because of this design, the algorithm is based purely on
genotype calls and heterozygote frequencies. Taking both of these into
account, the LOH Score algorithm is a log odds ratio of the probabilityof a region exhibiting LOH versus not exhibiting LOH. Because
levels of LOH can vary from sample to sample and region to region,
it is difficult to assign LOH score thresholds that always positively
identify regions exhibiting LOH. However, the LOH score is a valuablecalculation that can be used to detect chromosomal aberrations.

An odds ratio is defined as the ratio of the number of subjects in agroup with an event to the number of subjects without an event. The
log odds ratio for each of these hypotheses is computed using thecluster file to estimate heterozygote allele frequencies for every SNP,
assuming that genotyping calls are independent.
The LOH Score algorithm does not incorporate haplotype structures
and assumes that heterozygote frequencies in the training set are
representative of frequencies in the population under study. Therefore,
the LOH score is a generalization of what may be occurring in a regionof interest. In single-sample mode, the reference is a cluster file and it
is possible that copy-neutral LOH detected may be due to haplotype
block structure in the data.
A diverse panel of 270 HapMap samples, including Caucasian,
Han Chinese, Japanese, and Yoruba HapMap populations, isused to create the default cluster file and to calculate heterozygote
frequencies. Heterozygosity rates are estimated for the combinedgroup. If the population under study is not represented well by thisgroup, it is beneficial to create a new cluster file based on the datafrom the unique population. Independent of the platform used, false
positives in the LOH score may be due to some SNPs being rare in thestudied population but common in the diversity panel used to create
the cluster file.

LOH Score Example
To illustrate the flow of this algorithm, consider this simplified example.
For a window containing N SNPs, all resulting in homozygote callswith heterozygote frequencies of h = 0.1, let the genotyping error be
e = 0.001. The likelihood of LOH occurring is (1-e)N. The likelihood of no
LOH is (1-h)N. Therefore, the log odds ratio is: log10{(1-e)N/(1-h)N}, whichis the same as N{log10(1-e) -log10(1-h)}, which equals {N(h-e)/2.3}.
The odds of LOH or No LOH grows in a roughly linear fashion with thenumber of consecutive homozygotes.
Conversely, if that stretch contains M heterozygote calls, the likelihood of
LOH decreases and the equation is adjusted to (1-e)(N-M) × eM, because
heterozygotes in a region with LOH occur only through genotyping
errors. The likelihood of No LOH also changes and becomes(1-h)(N-M) x hM. The log odds ratio now becomes equal to{(N-M)(h-e)/2.3} + M{log(e)-log(h)}, which equals {(N-M)(h-e)/2.3}-2M.
In this case, the odds have diminished. When roughly 1 in 10 SNPs
receives a heterozygote call, the odds of both hypotheses are equal.
If more heterozygote calls are produced, the log odds ratio becomes
negative as it becomes less likely that these observations come from a
region with LOH.
It is important to remember that there are usually unknown haplotype
structures and population-dependent heterozygote frequencies thatmay play a role in the accuracy of the LOH score. However, this scoreis provided as a starting point to determine whether a particular stretchof homozygotes contains LOH.

Technical Note: Informatics

Summary
GenomeStudio provides several methods to analyze SNP andprobe intensity data to identify chromosomal regions with LOH
and copy number variations. The software plug-ins described inthis technical note are freely available to GenomeStudio users toprovide extended functionality.
Automated bookmarking algorithms save time by automaticallyscanning and categorizing samples. Researchers can use cnvPartition
to find and calculate copy numbers, or Homozygosity Detector to
identify extended tracts of LOH.
The LOH Score algorithm provides statistical information aboutchromosomal aberrations of interest. This information includes theprobability of LOH existing. This algorithm can be used to identify
interesting regions in large sample sets quickly or to analyze a morerefined region further.
The open architecture of Illumina GenomeStudio software allows forcustomized and advanced analysis tools for the downstream analysis
of Illumina DNA Analysis BeadChip Genotyping data. The plug-ins
described in this document can be downloaded from the Illumina
support page: http://support.illumina.com/downloads/genomestudio2-0-plug-ins.html.

References1.

pngu.mgh.harvard.edu/~purcell/plink/2.

Olshen AB, Venkatraman ES, Lucito R, Wigler M (2004) Circular binary
segmentation for the analysis of array-based DNA copy number data.
Biostatistics 5: 557–572.

3.

Venkatraman ES, Olshen AB (2007) A faster circular binary segmentation
algorithm for the analysis of array CGH data. Bioinformatics 23: 657–663.

Technical Note: Informatics

Illumina • 1.800.809.4566 toll-free (U.S.) • +1.858.202.4566 tel • techsupport@illumina.com • www.illumina.com
FOR RESEARCH USE ONLY
© 2007–2014 Illumina, Inc. All rights reserved.
Illumina, BeadArray, GenomeStudio, Infinium, the pumpkin orange color, and the Genetic Energy streaming bases design
are trademarks or registered trademarks of Illumina, Inc. All other brands and names contained herein are the propertyof their respective owners.
Pub. No. 970-2007-008 Current as of 04 January 2017