Skip to content

Commit d4c5a00

Browse files
committed
Add combine_c_region_dbs() to combine two existing C-region dbs
1 parent 156fdce commit d4c5a00

12 files changed

Lines changed: 216 additions & 94 deletions

‎NAMESPACE‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -184,6 +184,7 @@ export(
184184

185185
## combine_germline_dbs.R:
186186
combine_germline_dbs,
187+
combine_c_region_dbs,
187188

188189
## augment_germline_db.R:
189190
augment_germline_db_V,

‎R/combine_germline_dbs.R‎

Lines changed: 123 additions & 39 deletions
Original file line numberDiff line numberDiff line change
@@ -1,27 +1,65 @@
11
### =========================================================================
2-
### combine_germline_dbs()
2+
### combine_germline_dbs() & combine_c_region_dbs()
33
### -------------------------------------------------------------------------
44
###
55

66

7+
.check_comb_db_name <- function(db_name)
8+
{
9+
if (!isSingleNonWhiteString(db_name))
10+
stop(wmsg("'db_name' must be a single (non-empty) string"))
11+
if (has_whitespace(db_name))
12+
stop(wmsg("'db_name' cannot contain whitespace characters"))
13+
if (!has_prefix(db_name, "comb"))
14+
stop(wmsg("'db_name' must start with \"comb\""))
15+
}
16+
17+
.check_suffixes <- function(suffix1, suffix2)
18+
{
19+
if (!isSingleString(suffix1))
20+
stop(wmsg("'suffix1' must be a single string"))
21+
if (has_whitespace(suffix1))
22+
stop(wmsg("'suffix1' cannot contain whitespace characters"))
23+
if (!isSingleString(suffix2))
24+
stop(wmsg("'suffix2' must be a single string"))
25+
if (has_whitespace(suffix2))
26+
stop(wmsg("'suffix2' cannot contain whitespace characters"))
27+
if (suffix1 == suffix2)
28+
stop(wmsg("'suffix1' and 'suffix2' must be different"))
29+
}
30+
31+
.stop_on_colliding_allele_names <- function()
32+
{
33+
msg1 <- "Allele name collisions!"
34+
msg2 <- c("Even after adding 'suffix1' to the allele names of the first ",
35+
"db, and 'suffix2' to the allele names of the second db, some ",
36+
"allele names are still colliding. This can be avoided by using ",
37+
"more specific suffixes.")
38+
stop(wmsg(msg1), "\n ", wmsg(msg2))
39+
}
40+
741
.combine_locus_alleles <- function(locus, alleles1, alleles2, suffix1, suffix2)
842
{
943
stopifnot(isSingleNonWhiteString(locus),
1044
is(alleles1, "DNAStringSet"), is(alleles2, "DNAStringSet"),
11-
isSingleNonWhiteString(suffix1), isSingleNonWhiteString(suffix2))
45+
isSingleString(suffix1), isSingleString(suffix2))
1246
locus_alleles1 <- alleles1[substr(names(alleles1), 1L, 3L) == locus]
1347
locus_alleles2 <- alleles2[substr(names(alleles2), 1L, 3L) == locus]
1448
## Add user-supplied suffixes to disambiguate allele names.
1549
names(locus_alleles1) <- add_suffix(names(locus_alleles1), suffix1)
1650
names(locus_alleles2) <- add_suffix(names(locus_alleles2), suffix2)
17-
c(locus_alleles1, locus_alleles2)
51+
locus_alleles <- c(locus_alleles1, locus_alleles2)
52+
if (anyDuplicated(names(locus_alleles)))
53+
.stop_on_colliding_allele_names()
54+
locus_alleles
1855
}
1956

20-
.combine_region_dbs <- function(destdir, db_name1, db_name2, suffix1, suffix2,
57+
### We only support combining region dbs for the IG loci at the moment.
58+
.combine_region_dbs <- function(destdir, db_path1, db_path2, suffix1, suffix2,
2159
region_type, verbose=FALSE)
2260
{
23-
alleles1 <- load_germline_db(db_name1, region_types=region_type)
24-
alleles2 <- load_germline_db(db_name2, region_types=region_type)
61+
alleles1 <- readDNAStringSet(get_db_fasta_file(db_path1, region_type))
62+
alleles2 <- readDNAStringSet(get_db_fasta_file(db_path2, region_type))
2563
loci1 <- substr(names(alleles1), 1, 3)
2664
loci2 <- substr(names(alleles2), 1, 3)
2765
stopifnot(all(loci1 %in% IG_LOCI), all(loci2 %in% IG_LOCI))
@@ -33,21 +71,36 @@
3371

3472
locus2alleles <- locus2alleles[lengths(locus2alleles) != 0L]
3573
stopifnot(length(locus2alleles) != 0L)
36-
names(locus2alleles) <- paste0(names(locus2alleles), region_type, ".fasta")
3774

38-
fasta_files <- vapply(locus2alleles,
39-
function(alleles) {
75+
loci <- names(locus2alleles)
76+
fasta_files <- vapply(loci,
77+
function(locus) {
78+
alleles <- locus2alleles[[locus]]
4079
fasta_file <- tempfile(fileext=".fasta")
4180
writeXStringSet(alleles, fasta_file)
81+
if (verbose) {
82+
msg <- c("Combined and saved ", length(alleles), " ",
83+
region_type, " alleles for locus ", locus, " to ",
84+
"temp file ", basename(fasta_file))
85+
message(wmsg(msg))
86+
}
4287
fasta_file
4388
},
4489
character(1))
90+
if (verbose)
91+
message("")
92+
on.exit(unlink(fasta_files))
4593

94+
names(fasta_files) <- paste0(loci, region_type, ".fasta")
4695
create_region_db(fasta_files, destdir, region_type=region_type,
4796
verbose=verbose)
48-
lapply(fasta_files, nuke_file)
4997
}
5098

99+
100+
### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
101+
### combine_germline_dbs()
102+
###
103+
51104
.combine_db_intdata <- function(destdir, db_name1, db_name2, suffix1, suffix2,
52105
verbose=FALSE)
53106
{
@@ -81,49 +134,32 @@
81134
}
82135

83136
.do_combine_germline_dbs <- function(destdir,
84-
db_name1, db_name2, suffix1, suffix2,
137+
db_path1, db_path2, suffix1, suffix2,
85138
verbose=FALSE)
86139
{
87-
if (!isSingleNonWhiteString(db_name1))
88-
stop(wmsg("'db_name1' must be a single (non-empty) string"))
89-
if (!isSingleNonWhiteString(db_name2))
90-
stop(wmsg("'db_name2' must be a single (non-empty) string"))
91-
if (!isSingleString(suffix1))
92-
stop(wmsg("'suffix1' must be a single string"))
93-
if (has_whitespace(suffix1))
94-
stop(wmsg("'suffix1' cannot contain whitespace characters"))
95-
if (!isSingleString(suffix2))
96-
stop(wmsg("'suffix2' must be a single string"))
97-
if (has_whitespace(suffix2))
98-
stop(wmsg("'suffix2' cannot contain whitespace characters"))
99140
for (region_type in VDJ_REGION_TYPES)
100-
.combine_region_dbs(destdir, db_name1, db_name2, suffix1, suffix2,
141+
.combine_region_dbs(destdir, db_path1, db_path2, suffix1, suffix2,
101142
region_type, verbose=verbose)
143+
db_name1 <- basename(db_path1)
144+
db_name2 <- basename(db_path2)
102145
.combine_db_intdata(destdir, db_name1, db_name2, suffix1, suffix2,
103146
verbose=verbose)
104147
.combine_db_auxdata(destdir, db_name1, db_name2, suffix1, suffix2,
105148
verbose=verbose)
106149
}
107150

108-
109-
### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
110-
### combine_germline_dbs()
111-
###
112-
113-
### Combine two existing germline dbs, typically from two different
114-
### organisms, into one hybrid germline db.
151+
### Combine two existing germline dbs, typically (but not necessarily) from
152+
### two different organisms, into one hybrid germline db.
115153
### Note that allele names are usually not unique across species.
116154
### Arguments 'suffix1' and 'suffix2' can be used to disambiguate them.
117155
combine_germline_dbs <- function(db_name,
118156
db_name1, db_name2, suffix1, suffix2,
119157
overwrite=FALSE, verbose=FALSE)
120158
{
121-
if (!isSingleNonWhiteString(db_name))
122-
stop(wmsg("'db_name' must be a single (non-empty) string"))
123-
if (has_whitespace(db_name))
124-
stop(wmsg("'db_name' cannot contain whitespace characters"))
125-
if (!has_prefix(db_name, "comb"))
126-
stop(wmsg("'db_name' must start with \"comb\""))
159+
.check_comb_db_name(db_name)
160+
check_germline_db_name(db_name1, what="'db_name1'")
161+
check_germline_db_name(db_name2, what="'db_name2'")
162+
.check_suffixes(suffix1, suffix2)
127163

128164
if (!isTRUEorFALSE(overwrite))
129165
stop(wmsg("'overwrite' must be TRUE or FALSE"))
@@ -132,7 +168,10 @@ combine_germline_dbs <- function(db_name,
132168

133169
db_path <- get_germline_db_path(db_name)
134170
if (dir.exists(db_path) && !overwrite)
135-
stop_on_existing_cached_germline_db(db_name)
171+
stop_on_existing_germline_db(db_name)
172+
173+
db_path1 <- get_germline_db_path(db_name1)
174+
db_path2 <- get_germline_db_path(db_name2)
136175

137176
## We first create the combined germline db in a temporary folder, and,
138177
## only if successful, we replace 'db_path' with the temporary folder.
@@ -142,10 +181,55 @@ combine_germline_dbs <- function(db_name,
142181
tmp_destdir <- tempfile("germline_db_")
143182
dir.create(tmp_destdir)
144183
on.exit(nuke_file(tmp_destdir))
145-
.do_combine_germline_dbs(tmp_destdir, db_name1, db_name2, suffix1, suffix2,
184+
.do_combine_germline_dbs(tmp_destdir, db_path1, db_path2, suffix1, suffix2,
146185
verbose=verbose)
147186
rename_file(tmp_destdir, db_path, replace=TRUE)
148187

149188
invisible(db_name)
150189
}
151190

191+
192+
### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
193+
### combine_c_region_dbs()
194+
###
195+
196+
### Combine two existing C-region dbs, typically (but not necessarily) from
197+
### two different organisms, into one hybrid C-region db.
198+
### Note that allele names are usually not unique across species.
199+
### Arguments 'suffix1' and 'suffix2' can be used to disambiguate them.
200+
combine_c_region_dbs <- function(db_name,
201+
db_name1, db_name2, suffix1, suffix2,
202+
overwrite=FALSE, verbose=FALSE)
203+
{
204+
.check_comb_db_name(db_name)
205+
check_c_region_db_name(db_name1, what="'db_name1'")
206+
check_c_region_db_name(db_name2, what="'db_name2'")
207+
.check_suffixes(suffix1, suffix2)
208+
209+
if (!isTRUEorFALSE(overwrite))
210+
stop(wmsg("'overwrite' must be TRUE or FALSE"))
211+
if (!isTRUEorFALSE(verbose))
212+
stop(wmsg("'verbose' must be TRUE or FALSE"))
213+
214+
db_path <- get_c_region_db_path(db_name)
215+
if (dir.exists(db_path) && !overwrite)
216+
stop_on_existing_c_region_db(db_name)
217+
218+
db_path1 <- get_c_region_db_path(db_name1)
219+
db_path2 <- get_c_region_db_path(db_name2)
220+
221+
## We first create the combined C-region db in a temporary folder, and,
222+
## only if successful, we replace 'db_path' with the temporary folder.
223+
## Otherwise we destroy the temporary folder and raise an error. This
224+
## achieves atomicity and avoids loosing the content of the existing
225+
## 'db_path' in case something goes wrong.
226+
tmp_destdir <- tempfile("c_region_db_")
227+
dir.create(tmp_destdir)
228+
on.exit(nuke_file(tmp_destdir))
229+
.combine_region_dbs(tmp_destdir, db_path1, db_path2, suffix1, suffix2,
230+
"C", verbose=verbose)
231+
rename_file(tmp_destdir, db_path, replace=TRUE)
232+
233+
invisible(db_name)
234+
}
235+

‎R/create_c_region_db.R‎

Lines changed: 15 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -6,6 +6,20 @@
66
###
77

88

9+
.LIST_C_REGION_DB_TIP <- c(
10+
"Use list_c_region_dbs() to list all the C-region dbs ",
11+
"currently installed in the cache (see '?list_c_region_dbs')."
12+
)
13+
14+
stop_on_existing_c_region_db <- function(db_name)
15+
{
16+
msg1 <- c("C-region db ", db_name, " is already installed ",
17+
"in igblastr's persistent cache.")
18+
msg3 <- c("Use 'overwrite=TRUE' to reinstall.")
19+
stop(wmsg(msg1), "\n ", wmsg(.LIST_C_REGION_DB_TIP), "\n ", wmsg(msg3))
20+
}
21+
22+
923
### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1024
### list_loci_in_c_region_fasta_dir()
1125
###
@@ -59,16 +73,6 @@ list_loci_in_c_region_fasta_dir <- function(fasta_dir, loci_prefix)
5973
### create_c_region_db()
6074
###
6175

62-
.stop_on_existing_c_region_db <- function(destdir)
63-
{
64-
db_name <- basename(destdir)
65-
msg1 <- c("C-region db ", db_name, " is already installed.")
66-
msg2 <- c("Use list_c_region_dbs() to list the C-region databases ",
67-
"already installed on your machine (see '?list_c_region_dbs').")
68-
msg3 <- c("Use 'overwrite=TRUE' to reinstall.")
69-
stop(wmsg(msg1), "\n ", wmsg(msg2), "\n ", wmsg(msg3))
70-
}
71-
7276
### Creates a C-region db (constant regions) from a collection of FASTA
7377
### files (typically obtained from IMGT) for a given organism.
7478
### Note that 'destdir' will typically be the path to a subdir of the
@@ -93,7 +97,7 @@ create_c_region_db <- function(fasta_dir, loci, destdir,
9397
}
9498

9599
if (dir.exists(destdir) && !overwrite)
96-
.stop_on_existing_c_region_db(destdir)
100+
stop_on_existing_c_region_db(basename(destdir))
97101

98102
fasta_files <- .collect_C_fasta_files(fasta_dir, loci)
99103

‎R/create_region_db.R‎

Lines changed: 4 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -62,7 +62,7 @@
6262
.load_allele_set <- function(fasta_files, region_type,
6363
excluded_alleles=character(0), verbose=FALSE)
6464
{
65-
.checkarg_fasta_files(fasta_files)
65+
loci <- .infer_loci_from_filenames(fasta_files, region_type)
6666
stopifnot(is.character(excluded_alleles), isTRUEorFALSE(verbose))
6767

6868
if (verbose) {
@@ -71,7 +71,6 @@
7171
message(" o ", wmsg(msg, margin=4L), " ", appendLF=FALSE)
7272
}
7373

74-
loci <- .infer_loci_from_filenames(fasta_files, region_type)
7574
allele_sets <- lapply(seq_along(fasta_files),
7675
function(i) {
7776
allele_set <- readDNAStringSet(fasta_files[[i]])
@@ -284,7 +283,7 @@
284283
.write_region_db <-
285284
function(cleaned_allele_set, destdir, region_type,
286285
fasta_files, DORsummary, excluded_alleles=character(0),
287-
verbose=FALSE)
286+
verbose=FALSE)
288287
{
289288
stopifnot(is.data.frame(DORsummary))
290289
suffix <- DORsummary[ , "suffix"]
@@ -351,7 +350,7 @@ create_region_db <- function(fasta_files, destdir,
351350
disambiguate.allele.names=disambiguate.allele.names,
352351
summary.only=TRUE)
353352
.write_region_db(cleaned_allele_set, destdir, region_type,
354-
fasta_files, DORsummary, excluded_alleles=excluded_alleles,
353+
fasta_files, DORsummary, excluded_alleles=excluded_alleles,
355354
verbose=verbose)
356355
}
357356

@@ -412,7 +411,7 @@ create_V_region_db <- function(fasta_files, destdir,
412411
disambiguate.allele.names=disambiguate.allele.names,
413412
summary.only=TRUE)
414413
.write_region_db(cleaned_allele_set, destdir, "V",
415-
fasta_files, DORsummary,
414+
fasta_files, DORsummary,
416415
verbose=verbose)
417416
if (with.intdata) {
418417
if (verbose)

‎R/install_custom_germline_db.R‎

Lines changed: 7 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -3,23 +3,24 @@
33
### -------------------------------------------------------------------------
44

55

6-
### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
7-
### install_germline_db()
8-
###
9-
106
.LIST_GERMLINE_DB_TIP <- c(
117
"Use list_germline_dbs() to list all the germline dbs ",
128
"currently installed in the cache (see '?list_germline_dbs')."
139
)
1410

15-
stop_on_existing_cached_germline_db <- function(db_name)
11+
stop_on_existing_germline_db <- function(db_name)
1612
{
1713
msg1 <- c("Germline db ", db_name, " is already installed ",
1814
"in igblastr's persistent cache.")
1915
msg3 <- c("Use 'overwrite=TRUE' to reinstall.")
2016
stop(wmsg(msg1), "\n ", wmsg(.LIST_GERMLINE_DB_TIP), "\n ", wmsg(msg3))
2117
}
2218

19+
20+
### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
21+
### install_germline_db()
22+
###
23+
2324
.install_germline_db_succeeded <- function(db_name)
2425
{
2526
message(wmsg("Germline db ", db_name, " successfully installed ",
@@ -86,7 +87,7 @@ install_germline_db <- function(install_dir, db_name, fasta_dir, loci,
8687
if (if.exists == "no-op")
8788
return()
8889
if (if.exists == "error")
89-
stop_on_existing_cached_germline_db(db_name)
90+
stop_on_existing_germline_db(db_name)
9091
}
9192

9293
if (verbose) {

0 commit comments

Comments
 (0)