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
|
/*
* 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;
import beagleutil.ChromInterval;
import java.io.File;
/**
* <p>Interface {@code GeneticMap} represents a genetic map for one or more
* chromosomes.
* </p>
* <p>Instances of class {@code GeneticMap} are immutable.
* </p>
*
* @author Brian L. Browning {@code <browning@uw.edu>}
*/
public interface GeneticMap {
/**
* Returns the base position corresponding to the specified genetic map
* position. If the genetic position is not a map position then the base
* position is estimated from the nearest genetic map positions using
* linear interpolation.
*
* @param chrom the chromosome index
* @param geneticPosition the genetic position on the chromosome
* @return the base position corresponding to the specified genetic map
* position
* @throws IllegalArgumentException if the calculated base position
* exceeds {@code Integer.MAX_VALUE}
* @throws IllegalArgumentException if this genetic map has no
* map positions for the specified chromosome
* @throws IndexOutOfBoundsException if
* {@code chrom < 0 || chrom >= ChromIds.instance().size()}
*/
int basePos(int chrom, double geneticPosition);
/**
* Returns the genetic map position of the specified marker. The
* genetic map position is estimated using linear interpolation.
*
* @param marker a genetic marker
* @return the genetic map position of the specified marker
* @throws IllegalArgumentException if this genetic map has no
* map positions for the specified chromosome
* @throws NullPointerException if {@code marker == null}
*/
double genPos(Marker marker);
/**
* Returns the genetic map position of the specified genome coordinate.
* The genetic map position is estimated using linear interpolation.
*
* @param chrom the chromosome index
* @param basePosition the base coordinate on the chromosome
* @return the genetic map position of the specified genome coordinate
* @throws IllegalArgumentException if this genetic map has no
* map positions for the specified chromosome
* @throws IndexOutOfBoundsException if
* {@code chrom < 0 || chrom >= ChromIds.instance().size()}
*/
double genPos(int chrom, int basePosition);
/**
* Returns a string representation of this genetic map. The exact details
* of the representation are unspecified and subject to change.
*
* @return a string representation of this genetic map
*/
@Override
String toString();
/**
* Constructs and returns a genetic map from the specified data.
* If the specified map file is {@code null}, the returned genetic map
* will convert genome coordinates to genetic units by dividing by
* 1,000,000. If {@code (chromInt != null)} the genetic map will
* be restricted to chromosome {@code chromInt.chrom()}.
* @param file a PLINK-format genetic map file with cM units
* @param chromInt a chromosome interval
* @return a genetic map from the specified data.
* @throws IllegalArgumentException if any map position is infinite
* or {@code NaN}
* @throws NumberFormatException if the base position on any line of the map
* file is not a parsable integer
* @throws NumberFormatException if the genetic map position on any
* line of the map file is not a parsable double
* @throws IllegalArgumentException if a non-empty line of the specified
* genetic map file does not contain 4 fields
* @throws IllegalArgumentException if the map positions on each
* chromosome are not sorted in ascending order
* @throws IllegalArgumentException if there are duplicate
* base positions on a chromosome
* @throws IllegalArgumentException if all base positions on a chromosome
* have the same genetic map position
*/
static GeneticMap geneticMap(File file, ChromInterval chromInt) {
if (file==null) {
double scaleFactor = 1e-6;
return new PositionMap(scaleFactor);
}
else {
if (chromInt==null) {
return PlinkGenMap.fromPlinkMapFile(file);
}
else {
return PlinkGenMap.fromPlinkMapFile(file, chromInt.chrom());
}
}
}
/**
* Returns the an array of length {@code markers.nMarkers()} whose
* whose {@code j}-th element is the genetic map position of the
* {@code j}-th marker.
* @param genMap the genetic map
* @param markers the list of markers
* @return an array of genetic map positions
* @throws IllegalArgumentException if
* {@code markers.marker(0).chromIndex() != markers.marker(markers.nMarkers() - 1).chromIndex()}
* @throws IllegalArgumentException if the specified genetic map has no
* map positions for the specified chromosome
* @throws NullPointerException if {@code genMap == null || markers == null}
*/
static double[] genPos(GeneticMap genMap, Markers markers) {
if (markers.marker(0).chromIndex()
!= markers.marker(markers.size()-1).chromIndex()) {
throw new IllegalArgumentException("inconsistent data");
}
double[] genPos = new double[markers.size()];
for (int j=0; j<genPos.length; ++j) {
genPos[j] = genMap.genPos(markers.marker(j));
}
return genPos;
}
/**
* Returns the an array of length {@code markers.nMarkers()} whose
* whose {@code j}-th element is the genetic map position
* of the {@code j}-th marker.
* @param genMap the genetic map in cM units
* @param minGenDist the required minimum cM distance between successive
* markers
* @param markers the list of markers
* @return an array of genetic map positions
* @throws IllegalArgumentException if
* {@code markers.marker(0).chromIndex() != markers.marker(markers.nMarkers() - 1).chromIndex()}
* @throws IllegalArgumentException if the specified genetic map has no
* map positions for the specified chromosome
* @throws IllegalArgumentException if {@code Double.isFinite(minDist) == false}
* @throws NullPointerException if {@code genMap == null || markers == null}
*/
static double[] genPos(GeneticMap genMap, double minGenDist, Markers markers) {
if (markers.marker(0).chromIndex()
!= markers.marker(markers.size()-1).chromIndex()) {
throw new IllegalArgumentException("inconsistent data");
}
if (Double.isFinite(minGenDist)==false) {
throw new IllegalArgumentException(String.valueOf(minGenDist));
}
double[] genPos = new double[markers.size()];
genPos[0] = genMap.genPos(markers.marker(0));
double lastMapPos = genPos[0];
for (int j=1; j<genPos.length; ++j) {
double mapPos = genMap.genPos(markers.marker(j));
double dist = Math.max((mapPos - lastMapPos), minGenDist);
genPos[j] = genPos[j-1] + dist;
lastMapPos = mapPos;
}
return genPos;
}
}
|
