1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
|
/*
* Copyright (C) 2014-2021 Brian L. Browning
*
* This file is part of Beagle
*
* Beagle is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* Beagle is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*/
package vcf;
/**
* <p>Interface {@code GT} represents genotype data
* for a list of markers and a list of samples.
* </p>
* <p>All instances of {@code GT} are required to be immutable.
* </p>
*
* @author Brian L. Browning {@code <browning@uw.edu>}
*/
public interface GT {
/**
* Returns {@code true} if the markers are ordered by decreasing chromosome
* base position, and returns {@code false} otherwise.
* @return {@code true} if the markers are ordered by decreasing chromosome
* base position
*/
boolean isReversed();
/**
* Returns the number of markers.
* @return the number of markers
*/
int nMarkers();
/**
* Returns the specified marker.
* @param marker a marker index
* @return the specified marker
* @throws IndexOutOfBoundsException if
* {@code marker < 0 || marker >= this.nMarkers()}
*/
Marker marker(int marker);
/**
* Returns the list of markers in order of increasing chromosome position.
* If {@code (this.isReversed() == false)} then
* {@code (this.markers().marker(j).equals(this.marker(j)) == true)} for
* all {@code (0 <= j && j < this.nMarkers())}.
* If {@code (this.isReversed() == true)} then
* {@code (this.markers().marker(this.nMarkers() - 1 - j).equals(this.marker(j)) == true)}
* for all {@code (0 <= j && j < this.nMarkers())}
* @return the list of markers in order of increasing chromosome position
*/
Markers markers();
/**
* Returns the number of haplotypes. The returned value is equal to
* {@code 2*this.nSamples()}.
* @return the number of haplotypes
*/
int nHaps();
/**
* Returns the number of samples.
* @return the number of samples
*/
int nSamples();
/**
* Returns the list of samples.
* @return the list of samples
*/
Samples samples();
/**
* Returns {@code true} if the genotype for each marker and sample
* has non-missing alleles and is either haploid or diploid with
* a phased allele separator, and returns {@code false} otherwise.
* @return {@code true} if the genotype for each marker and sample
* is a phased, non-missing genotype
*/
boolean isPhased();
/**
* Returns the first allele for the specified marker and sample
* or return -1 if the allele is missing. The order of the two alleles
* is unspecified if {@code this.isPhased() == false}.
* @param marker the marker index
* @param sample the sample index
* @return the first allele for the specified marker and sample
*
* @throws IndexOutOfBoundsException if
* {@code marker < 0 || marker >= this.nMarkers()}
* @throws IndexOutOfBoundsException if
* {@code sample < 0 || sample >= this.nSamples()}
*/
int allele1(int marker, int sample);
/**
* Returns the second allele for the specified marker and sample
* or return -1 if the allele is missing. The order of the two alleles
* is unspecified if {@code this.isPhased() == false}.
* @param marker the marker index
* @param sample the sample index
* @return the allele for the specified marker and sample
*
* @throws IndexOutOfBoundsException if
* {@code marker < 0 || marker >= this.nMarkers()}
* @throws IndexOutOfBoundsException if
* {@code sample < 0 || sample >= this.nSamples()}
*/
int allele2(int marker, int sample);
/**
* Returns the allele on the specified haplotype for the specified marker
* or return -1 if the allele is missing. The order of the two alleles
* is unspecified if {@code this.isPhased() == false}.
* @param marker the marker index
* @param hap the haplotype index
* @return the allele on the specified haplotype for the specified marker
*
* @throws IndexOutOfBoundsException if
* {@code marker < 0 || marker >= this.nMarkers()}
* @throws IndexOutOfBoundsException if
* {@code hap < 0 || hap >= this.nHaps()}
*/
int allele(int marker, int hap);
/**
* Returns a {@code GT} instance restricted to genotype data for
* the specified markers.
* @param markers the list of markers in the returned instance
* @param indices a list of distinct marker indices (from
* {@code this.markers())} in increasing order
* @return a {@code GT} instance restricted to genotype data for
* the specified markers
*
* @throws IndexOutOfBoundsException if there exists {@code j} such that
* {@code (0 <= j && j < indices.length)} such that
* {@code (indices[j] < 0 || indices[j] >= this.nMarkers())}
* @throws IllegalArgumentException if there exists {@code j} such that
* {@code (1 <= j && j < indices.length)} such that
* {@code (indices[j] <= indice[j - 1])}
* @throws IllegalArgumentException if there exists {@code j} such that
* {@code (0 <= j && j < indices.length)} such that
* {@code (this.marker(indices[j]).equals(markers.marker(j)) == false)}
* @throws NullPointerException if {@code indices == null}
* @throws UnsupportedOperationException if {@code this.isReversed() == true}
*/
GT restrict(Markers markers, int[] indices);
/**
* Returns a new {@code GT} instance restricted to genotype data for
* the specified markers.
* @param start the start marker (inclusive)
* @param end the end marker (exclusive)
* @return a {@code GT} instance restricted to genotype data for
* the specified markers
* @throws IllegalArgumentException if {@code start >= end}
* @throws IndexOutOfBoundsException if
* {@code start < 0 || end > this.markers()}
*/
GT restrict(int start, int end);
/**
* Returns the probability of the observed data for the specified marker
* and sample if the specified pair of unordered alleles is the true
* genotype. Returns {@code 1.0f} if the corresponding genotype
* determined by the {@code isPhased()}, {@code allele1()}, and
* {@code allele2()} methods is consistent with the specified ordered
* genotype, and returns {@code 0.0f} otherwise.
*
* @param gt the genotype data
* @param marker the marker index
* @param sample the sample index
* @param allele1 the first allele index
* @param allele2 the second allele index
* @return the probability of the observed data for the specified marker
* and sample if the specified pair of ordered alleles is the true
* ordered genotype
*
* @throws IndexOutOfBoundsException if
* {@code marker < 0 || marker >= this.nMarkers()}
* @throws IndexOutOfBoundsException if
* {@code samples < 0 || samples >= this.nSamples()}
* @throws IndexOutOfBoundsException if
* {@code allele1 < 0 || allele1 >= this.marker(marker).nAlleles()}
* @throws IndexOutOfBoundsException if
* {@code allele2 < 0 || allele2 >= this.marker(marker).nAlleles()}
*/
static float gl(GT gt, int marker, int sample, int allele1, int allele2) {
int nAlleles = gt.marker(marker).nAlleles();
if (allele1 < 0 || allele1 >= nAlleles) {
String s = "invalid alleles: (" + allele1 + "): " + marker;
throw new IllegalArgumentException(s);
}
if (allele2 < 0 || allele2 >= nAlleles) {
String s = "invalid alleles: (" + allele2 + "): " + marker;
throw new IllegalArgumentException(s);
}
int a1 = gt.allele1(marker, sample);
int a2 = gt.allele2(marker, sample);
boolean consistent = (a1==-1 || a1==allele1) && (a2==-1 || a2==allele2);
if (consistent==false) {
consistent = (a1==-1 || a1==allele2) && (a2==-1 || a2==allele1);
}
return consistent ? 1.0f : 0.0f;
}
}
|
