Beginning with genotypes for the colony samples, we want to ready those data for self-assignment using rubias.
This workflow uses some functions from the CKMRsim package, which is available from https://github.com/eriqande/CKMRsim
Load data and libraries
library(tidyverse)
library(CKMRsim)
library(stringr)
meta <- readRDS("data/meta-data-tibble.rds")
genos <- readRDS("data/april_called_genos_na_explicit.rds") %>%
filter(NMFS_DNA_ID %in% meta$NMFS_DNA_ID) # drop those we don't have meta data for
samples <- readRDS("data/sample-sheet-tibble.rds") %>%
filter(NMFS_DNA_ID %in% meta$NMFS_DNA_ID)
Some initial filters
keepers <- meta %>%
dplyr::select(NMFS_DNA_ID) %>%
unlist() %>% unname()
Whittle the genotype data down to those from the meta data:
kgenos <- genos %>%
filter(NMFS_DNA_ID %in% keepers)
611 samples
Check regenotyped samples for consistency
There should be some of these from my re-genotyping effort with Chagulak.
# regeno <- kgenos %>%
# select(NMFS_DNA_ID) %>%
# group_by(NMFS_DNA_ID) %>%
# tally() %>%
# arrange(desc(n)) %>%
# filter(n > 500) %>% # these are the 61 that I re-genotyped.
# select(NMFS_DNA_ID)
Take highest read-depth call for multiply-genotyped DNA_IDs
I re-genotyped some samples with poor quality from the first sequencing run. Here, if an individual is multiply-genotyped, take the genotype with the highest total read depth.
# slow-ish function to get the total read depth column
tdepth <- function(a, d) {
if(any(is.na(a))) {
return(NA)
}
if(a[1]==a[2]) {
return(d[1])
} else {
return(d[1] + d[2])
}
}
# this takes the highest read-depth instance of each duplicately-genotyped individual.
geno_one_each <- kgenos %>%
group_by(NMFS_DNA_ID, locus, gtseq_run, id) %>%
mutate(total_depth = tdepth(allele, depth)) %>%
ungroup() %>%
arrange(NMFS_DNA_ID, locus, desc(total_depth), gtseq_run, id, depth) %>%
group_by(NMFS_DNA_ID, locus) %>%
mutate(rank = 1:n()) %>%
ungroup() %>%
filter(rank <= 2)
# geno_one_each %>%
# select(NMFS_DNA_ID) %>%
# unique()
Before tossing individuals with too much missing data, we want to implement our total read depth threshold.
# out of curiosity, how many samples would be affected by the read depth filter at this point?
how_many <- geno_one_each %>%
filter(total_depth < 20) %>% # only loci for which total read depth < 20
group_by(NMFS_DNA_ID) %>%
dplyr::select(NMFS_DNA_ID) %>%
count() %>%
mutate(n = n/2) %>% # to get to per-locus, rather than per-allele
#arrange(desc(n)) %>%
filter(n > 28)
how_many
There are 17 birds that will have missing data at > 28 loci because of the total depth threshold of 20 reads.
# if total depth >= 20, print the allele; if not, print NA
geno_one_each2 <- geno_one_each %>%
mutate(allele = ifelse(total_depth >=20, allele, NA))
Toss out indivs with more than 28 missing loci (20%)
Now, toss out any individual with more than 28 missing loci
no_hi_missers <- geno_one_each2 %>%
group_by(NMFS_DNA_ID) %>%
filter(sum(!is.na(allele)) >= (113*2))
# How many is that?
kept <- no_hi_missers %>%
dplyr::select(NMFS_DNA_ID) %>%
unique()
# save that to look at duplicated samples and other things in 03-CKMR-power-analysis.Rmd
no_hi_missers %>%
write_csv("csv_outputs/no_hi_missers.csv")
Making a GSI (Rubias) data set
For the baselines, we want to treat the colony locations as populations.
# rather than species, we want to hone in on the location information.
# Let's organize that a bit.
meta_select <- meta %>%
dplyr::select(NMFS_DNA_ID, BOX_ID, BOX_POSITION, SAMPLE_ID, BATCH_ID, PROJECT_NAME, GENUS, SPECIES, REPORTED_LIFE_STAGE, COLLECTION_DATE, `Marine::LOCATION_COMMENTS_M`)
# rename the location column
colnames(meta_select)[11] <- "LOCATION"
meta_select <- meta_select %>%
mutate(LOCATION = ifelse(LOCATION == "St. George", "StGeorge", LOCATION)) %>%
mutate(LOCATION = ifelse(LOCATION == "St. Paul", "StPaul", LOCATION)) %>%
mutate(LOCATION = ifelse(LOCATION == "St. Matthew", "StMatthew", LOCATION))
meta_keep <- meta_select %>%
group_by(LOCATION)
Assignment to main colonies
# modify the meta data colony identifiers to be consistent for the four major colonies
meta_keep <- meta_keep %>%
mutate(COLONY = ifelse(LOCATION == "StGeorge", "Pribilof",
ifelse(LOCATION == "StPaul", "Pribilof",
ifelse(LOCATION == "Chowiet", "Semidi",
ifelse(LOCATION == "Kateekuk", "Semidi", LOCATION)
)
)))
Let’s pick them all out:
nofu_ids <- meta_keep %>%
ungroup() %>%
mutate(group = COLONY) %>%
dplyr::select(NMFS_DNA_ID, group) %>%
arrange(group)
Allele frequencies
With those we can just filter down the genos to the ones that we want, and then we can get it into the format required for CKMR.
# we will use this some more
kg2 <- no_hi_missers %>%
select(NMFS_DNA_ID, locus, allele) %>%
mutate(Chrom = "GTseq") %>%
mutate(Pos = as.integer(factor(locus, levels = unique(locus)))) %>%
rename(Locus = locus,
Allele = allele) %>%
select(NMFS_DNA_ID, Chrom, Locus, Pos, Allele) %>%
ungroup()
# get the allele freqs
kg_ckmr_markers <- kg2 %>%
filter(!is.na(Allele)) %>% # it is vital to filter out the NAs at this stage
group_by(Chrom, Locus, Pos, Allele) %>%
summarise(counts = n()) %>%
group_by(Locus, Pos) %>%
mutate(Freq = counts / sum(counts)) %>%
select(-counts) %>%
mutate(AlleIdx = 1,
LocIdx = 1) %>%
reindex_markers(.)
Making genotype matrices & checking for duplicated samples
nofu_haps <- no_hi_missers %>%
filter(!is.na(allele)) %>% # once again, it is critical to remove these at this point
select(NMFS_DNA_ID, locus, gene_copy, allele) %>%
rename(Locus = locus, Allele = allele)
nofu_idx_frame <- kg_ckmr_markers %>%
select(Locus, Allele, LocIdx, AlleIdx) %>%
group_by(Locus) %>%
mutate(NumA = n()) %>% # get the number of alleles at each locus
ungroup() %>%
left_join(nofu_haps, .) %>% # join the alle_idx's onto the actual genotype data
select(NMFS_DNA_ID, Locus, gene_copy, LocIdx, NumA, AlleIdx) %>%
spread(key = gene_copy, value = AlleIdx) %>%
mutate(GenoIdx = index_ab(a = `1`, b = `2`, A = NumA))
Joining, by = c("Locus", "Allele")
# make a matrix of genotype integers
wide_nofu <- nofu_idx_frame %>%
select(NMFS_DNA_ID, LocIdx, GenoIdx) %>%
spread(data = ., key = LocIdx, value = GenoIdx)
Don’t forget to set NA’s to 0, and then decrease each value by 1:
nofu_mat <- as.matrix(wide_nofu[, -1])
rownames(nofu_mat) <- wide_nofu$NMFS_DNA_ID
nofu_mat[is.na(nofu_mat)] <- 0
nofu_mat <- nofu_mat - 1
storage.mode(nofu_mat) <- "integer"
Looking for duplicated samples
Looking back at Andy Ramey’s USGS meta data, there are some duplicate samples? Or at least birds that may have been sampled multiple times.
We can quickly look through rocky_mat for pairs of indivs with lots of matching genotypes.
matchers <- pairwise_geno_id(S = nofu_mat, max_miss = 12) %>%
arrange(num_mismatch) %>%
mutate(NMFS_DNA_ID_1 = rownames(nofu_mat)[ind1],
NMFS_DNA_ID_2 = rownames(nofu_mat)[ind2])
matchers %>%
arrange(desc(num_loc)) # need to use the filtered genotype data for this...
There are 13 samples that match perfectly. I should get rid of one of each pair.
Now, to deal with these duplicates, here is what I will do: since all of the matchers have 0 mismatched loci, we are going to take just one from each pairs.
def_same <- matchers %>%
filter(num_mismatch <= 3)
toss_these <- c(intersect(def_same$ind1, def_same$ind2), def_same$ind2) %>%
unique()
nofu_mat_tossed <- nofu_mat[-toss_these,]
no_hi_missers_clean <- no_hi_missers %>%
filter(NMFS_DNA_ID %in% rownames(nofu_mat_tossed))
# save that for the bycatch assignments
write_csv(no_hi_missers_clean, "csv_outputs/baseline_no_hi_missers_clean.csv")
LS0tCnRpdGxlOiAiMDEtdGlkeS1jb2xvbnktZ2Vub3R5cGVzIgpvdXRwdXQ6IGh0bWxfbm90ZWJvb2sKLS0tCgpCZWdpbm5pbmcgd2l0aCBnZW5vdHlwZXMgZm9yIHRoZSBjb2xvbnkgc2FtcGxlcywgd2Ugd2FudCB0byByZWFkeSB0aG9zZSBkYXRhIGZvciBzZWxmLWFzc2lnbm1lbnQgdXNpbmcgcnViaWFzLgoKVGhpcyB3b3JrZmxvdyB1c2VzIHNvbWUgZnVuY3Rpb25zIGZyb20gdGhlIENLTVJzaW0gcGFja2FnZSwgd2hpY2ggaXMgYXZhaWxhYmxlIGZyb20gaHR0cHM6Ly9naXRodWIuY29tL2VyaXFhbmRlL0NLTVJzaW0KCgoKTG9hZCBkYXRhIGFuZCBsaWJyYXJpZXMKYGBge3IgbG9hZC1kYXRhLWFuZC1saWJzfQpsaWJyYXJ5KHRpZHl2ZXJzZSkKbGlicmFyeShDS01Sc2ltKQpsaWJyYXJ5KHN0cmluZ3IpCgptZXRhIDwtIHJlYWRSRFMoImRhdGEvbWV0YS1kYXRhLXRpYmJsZS5yZHMiKQpnZW5vcyA8LSByZWFkUkRTKCJkYXRhL2FwcmlsX2NhbGxlZF9nZW5vc19uYV9leHBsaWNpdC5yZHMiKSAlPiUKICBmaWx0ZXIoTk1GU19ETkFfSUQgJWluJSBtZXRhJE5NRlNfRE5BX0lEKSAgIyBkcm9wIHRob3NlIHdlIGRvbid0IGhhdmUgbWV0YSBkYXRhIGZvcgpzYW1wbGVzIDwtIHJlYWRSRFMoImRhdGEvc2FtcGxlLXNoZWV0LXRpYmJsZS5yZHMiKSAlPiUKICBmaWx0ZXIoTk1GU19ETkFfSUQgJWluJSBtZXRhJE5NRlNfRE5BX0lEKQpgYGAKCgoKIyMgU29tZSBpbml0aWFsIGZpbHRlcnMKCmBgYHtyIGdldC1mdWxtYXJzfQprZWVwZXJzIDwtIG1ldGEgJT4lCiAgZHBseXI6OnNlbGVjdChOTUZTX0ROQV9JRCkgJT4lCiAgdW5saXN0KCkgJT4lIHVubmFtZSgpCiAgCmBgYAoKV2hpdHRsZSB0aGUgZ2Vub3R5cGUgZGF0YSBkb3duIHRvIHRob3NlIGZyb20gdGhlIG1ldGEgZGF0YToKYGBge3Igd2hpdHRsZS1nZW5vc30Ka2dlbm9zIDwtIGdlbm9zICU+JQogIGZpbHRlcihOTUZTX0ROQV9JRCAlaW4lIGtlZXBlcnMpCgpgYGAKNjExIHNhbXBsZXMKCiMjIENoZWNrIHJlZ2Vub3R5cGVkIHNhbXBsZXMgZm9yIGNvbnNpc3RlbmN5CgpUaGVyZSBzaG91bGQgYmUgc29tZSBvZiB0aGVzZSBmcm9tIG15IHJlLWdlbm90eXBpbmcgZWZmb3J0IHdpdGggQ2hhZ3VsYWsuCmBgYHtyfQojIHJlZ2VubyA8LSBrZ2Vub3MgJT4lCiMgICBzZWxlY3QoTk1GU19ETkFfSUQpICU+JQojICAgZ3JvdXBfYnkoTk1GU19ETkFfSUQpICU+JQojICAgdGFsbHkoKSAlPiUKIyAgIGFycmFuZ2UoZGVzYyhuKSkgJT4lCiMgICBmaWx0ZXIobiA+IDUwMCkgJT4lICMgdGhlc2UgYXJlIHRoZSA2MSB0aGF0IEkgcmUtZ2Vub3R5cGVkLgojICAgc2VsZWN0KE5NRlNfRE5BX0lEKQpgYGAKCgojIyMgVGFrZSBoaWdoZXN0IHJlYWQtZGVwdGggY2FsbCBmb3IgbXVsdGlwbHktZ2Vub3R5cGVkIEROQV9JRHMKCkkgcmUtZ2Vub3R5cGVkIHNvbWUgc2FtcGxlcyB3aXRoIHBvb3IgcXVhbGl0eSBmcm9tIHRoZSBmaXJzdCBzZXF1ZW5jaW5nIHJ1bi4gSGVyZSwgaWYgYW4gaW5kaXZpZHVhbCBpcyBtdWx0aXBseS1nZW5vdHlwZWQsIHRha2UgdGhlIGdlbm90eXBlIHdpdGggdGhlIGhpZ2hlc3QgdG90YWwgcmVhZCBkZXB0aC4gCmBgYHtyIHRha2UtanVzdC1vbmV9CiMgc2xvdy1pc2ggZnVuY3Rpb24gdG8gZ2V0IHRoZSB0b3RhbCByZWFkIGRlcHRoIGNvbHVtbgp0ZGVwdGggPC0gZnVuY3Rpb24oYSwgZCkgewogIGlmKGFueShpcy5uYShhKSkpIHsKICAgIHJldHVybihOQSkKICB9CiAgaWYoYVsxXT09YVsyXSkgewogICAgcmV0dXJuKGRbMV0pCiAgfSBlbHNlIHsKICAgIHJldHVybihkWzFdICsgZFsyXSkKICB9CiAgCn0KIyB0aGlzIHRha2VzIHRoZSBoaWdoZXN0IHJlYWQtZGVwdGggaW5zdGFuY2Ugb2YgZWFjaCBkdXBsaWNhdGVseS1nZW5vdHlwZWQgaW5kaXZpZHVhbC4KZ2Vub19vbmVfZWFjaCA8LSBrZ2Vub3MgJT4lCiAgZ3JvdXBfYnkoTk1GU19ETkFfSUQsIGxvY3VzLCBndHNlcV9ydW4sIGlkKSAlPiUKICBtdXRhdGUodG90YWxfZGVwdGggPSB0ZGVwdGgoYWxsZWxlLCBkZXB0aCkpICU+JQogIHVuZ3JvdXAoKSAlPiUKICBhcnJhbmdlKE5NRlNfRE5BX0lELCBsb2N1cywgZGVzYyh0b3RhbF9kZXB0aCksIGd0c2VxX3J1biwgaWQsIGRlcHRoKSAlPiUKICBncm91cF9ieShOTUZTX0ROQV9JRCwgbG9jdXMpICU+JQogIG11dGF0ZShyYW5rID0gMTpuKCkpICU+JQogIHVuZ3JvdXAoKSAlPiUKICBmaWx0ZXIocmFuayA8PSAyKQoKIyBnZW5vX29uZV9lYWNoICU+JQojICAgc2VsZWN0KE5NRlNfRE5BX0lEKSAlPiUKIyAgIHVuaXF1ZSgpCmBgYAoKCkJlZm9yZSB0b3NzaW5nIGluZGl2aWR1YWxzIHdpdGggdG9vIG11Y2ggbWlzc2luZyBkYXRhLCB3ZSB3YW50IHRvIGltcGxlbWVudCBvdXIgdG90YWwgcmVhZCBkZXB0aCB0aHJlc2hvbGQuCmBgYHtyIHRvdGFsLWRlcHRofQojIG91dCBvZiBjdXJpb3NpdHksIGhvdyBtYW55IHNhbXBsZXMgd291bGQgYmUgYWZmZWN0ZWQgYnkgdGhlIHJlYWQgZGVwdGggZmlsdGVyIGF0IHRoaXMgcG9pbnQ/Cmhvd19tYW55IDwtIGdlbm9fb25lX2VhY2ggJT4lCiAgZmlsdGVyKHRvdGFsX2RlcHRoIDwgMjApICU+JSAjIG9ubHkgbG9jaSBmb3Igd2hpY2ggdG90YWwgcmVhZCBkZXB0aCA8IDIwCiAgZ3JvdXBfYnkoTk1GU19ETkFfSUQpICU+JQogIGRwbHlyOjpzZWxlY3QoTk1GU19ETkFfSUQpICU+JQogIGNvdW50KCkgJT4lCiAgbXV0YXRlKG4gPSBuLzIpICU+JSAjIHRvIGdldCB0byBwZXItbG9jdXMsIHJhdGhlciB0aGFuIHBlci1hbGxlbGUKICAjYXJyYW5nZShkZXNjKG4pKSAlPiUKICBmaWx0ZXIobiA+IDI4KQoKaG93X21hbnkKYGBgClRoZXJlIGFyZSAxNyBiaXJkcyB0aGF0IHdpbGwgaGF2ZSBtaXNzaW5nIGRhdGEgYXQgPiAyOCBsb2NpIGJlY2F1c2Ugb2YgdGhlIHRvdGFsIGRlcHRoIHRocmVzaG9sZCBvZiAyMCByZWFkcy4KCgpgYGB7ciBmaWx0ZXItbG93LWRlcHRofQojIGlmIHRvdGFsIGRlcHRoID49IDIwLCBwcmludCB0aGUgYWxsZWxlOyBpZiBub3QsIHByaW50IE5BCmdlbm9fb25lX2VhY2gyIDwtIGdlbm9fb25lX2VhY2ggJT4lCiAgIG11dGF0ZShhbGxlbGUgPSBpZmVsc2UodG90YWxfZGVwdGggPj0yMCwgYWxsZWxlLCBOQSkpCgpgYGAKCgojIyMgVG9zcyBvdXQgaW5kaXZzIHdpdGggbW9yZSB0aGFuIDI4IG1pc3NpbmcgbG9jaSAoMjAlKQpOb3csIHRvc3Mgb3V0IGFueSBpbmRpdmlkdWFsIHdpdGggbW9yZSB0aGFuIDI4IG1pc3NpbmcgbG9jaQpgYGB7ciB0b3NzLW1pc3NlcnN9Cm5vX2hpX21pc3NlcnMgPC0gZ2Vub19vbmVfZWFjaDIgJT4lIAogIGdyb3VwX2J5KE5NRlNfRE5BX0lEKSAlPiUKICBmaWx0ZXIoc3VtKCFpcy5uYShhbGxlbGUpKSA+PSAoMTEzKjIpKQoKIyBIb3cgbWFueSBpcyB0aGF0PwprZXB0IDwtIG5vX2hpX21pc3NlcnMgJT4lCiAgZHBseXI6OnNlbGVjdChOTUZTX0ROQV9JRCkgJT4lCiAgdW5pcXVlKCkKCiMgc2F2ZSB0aGF0IHRvIGxvb2sgYXQgZHVwbGljYXRlZCBzYW1wbGVzIGFuZCBvdGhlciB0aGluZ3MgaW4gMDMtQ0tNUi1wb3dlci1hbmFseXNpcy5SbWQKbm9faGlfbWlzc2VycyAlPiUKICB3cml0ZV9jc3YoImNzdl9vdXRwdXRzL25vX2hpX21pc3NlcnMuY3N2IikKCmBgYAoKIyMgTWFraW5nIGEgR1NJIChSdWJpYXMpIGRhdGEgc2V0CgpGb3IgdGhlIGJhc2VsaW5lcywgd2Ugd2FudCB0byB0cmVhdCB0aGUgY29sb255IGxvY2F0aW9ucyBhcyBwb3B1bGF0aW9ucy4KYGBge3IgbWFrZS1ncm91cC1pbnRzfQojIHJhdGhlciB0aGFuIHNwZWNpZXMsIHdlIHdhbnQgdG8gaG9uZSBpbiBvbiB0aGUgbG9jYXRpb24gaW5mb3JtYXRpb24uIAojIExldCdzIG9yZ2FuaXplIHRoYXQgYSBiaXQuCm1ldGFfc2VsZWN0IDwtIG1ldGEgJT4lCiAgZHBseXI6OnNlbGVjdChOTUZTX0ROQV9JRCwgQk9YX0lELCBCT1hfUE9TSVRJT04sIFNBTVBMRV9JRCwgQkFUQ0hfSUQsIFBST0pFQ1RfTkFNRSwgR0VOVVMsIFNQRUNJRVMsIFJFUE9SVEVEX0xJRkVfU1RBR0UsIENPTExFQ1RJT05fREFURSwgYE1hcmluZTo6TE9DQVRJT05fQ09NTUVOVFNfTWApCgojIHJlbmFtZSB0aGUgbG9jYXRpb24gY29sdW1uCmNvbG5hbWVzKG1ldGFfc2VsZWN0KVsxMV0gPC0gIkxPQ0FUSU9OIgoKbWV0YV9zZWxlY3QgPC0gbWV0YV9zZWxlY3QgJT4lCiAgbXV0YXRlKExPQ0FUSU9OID0gaWZlbHNlKExPQ0FUSU9OID09ICJTdC4gR2VvcmdlIiwgIlN0R2VvcmdlIiwgTE9DQVRJT04pKSAlPiUKICBtdXRhdGUoTE9DQVRJT04gPSBpZmVsc2UoTE9DQVRJT04gPT0gIlN0LiBQYXVsIiwgIlN0UGF1bCIsIExPQ0FUSU9OKSkgJT4lCiAgbXV0YXRlKExPQ0FUSU9OID0gaWZlbHNlKExPQ0FUSU9OID09ICJTdC4gTWF0dGhldyIsICJTdE1hdHRoZXciLCBMT0NBVElPTikpIAoKbWV0YV9rZWVwIDwtIG1ldGFfc2VsZWN0ICU+JQogIGdyb3VwX2J5KExPQ0FUSU9OKSAKCmBgYAoKIyMgQXNzaWdubWVudCB0byBtYWluIGNvbG9uaWVzCgpgYGB7ciBtdXRhdGUtbWV0YX0KIyBtb2RpZnkgdGhlIG1ldGEgZGF0YSBjb2xvbnkgaWRlbnRpZmllcnMgdG8gYmUgY29uc2lzdGVudCBmb3IgdGhlIGZvdXIgbWFqb3IgY29sb25pZXMKbWV0YV9rZWVwIDwtIG1ldGFfa2VlcCAlPiUKICBtdXRhdGUoQ09MT05ZID0gaWZlbHNlKExPQ0FUSU9OID09ICJTdEdlb3JnZSIsICJQcmliaWxvZiIsCiAgICAgICAgICAgICAgICAgICAgICAgICBpZmVsc2UoTE9DQVRJT04gPT0gIlN0UGF1bCIsICJQcmliaWxvZiIsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgaWZlbHNlKExPQ0FUSU9OID09ICJDaG93aWV0IiwgIlNlbWlkaSIsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGlmZWxzZShMT0NBVElPTiA9PSAiS2F0ZWVrdWsiLCAiU2VtaWRpIiwgTE9DQVRJT04pCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICkKICAgICAgICAgICAgICAgICAgICAgICAgICkpKQpgYGAKCgpMZXQncyBwaWNrIHRoZW0gYWxsIG91dDoKYGBge3IgcGljay1vdXR9Cm5vZnVfaWRzIDwtIG1ldGFfa2VlcCAlPiUKICB1bmdyb3VwKCkgJT4lCiAgbXV0YXRlKGdyb3VwID0gQ09MT05ZKSAlPiUKICBkcGx5cjo6c2VsZWN0KE5NRlNfRE5BX0lELCBncm91cCkgJT4lCiAgYXJyYW5nZShncm91cCkKCmBgYAoKCiMjIEFsbGVsZSBmcmVxdWVuY2llcwpXaXRoIHRob3NlIHdlIGNhbiBqdXN0IGZpbHRlciBkb3duIHRoZSBnZW5vcyB0byB0aGUgb25lcyB0aGF0IHdlIHdhbnQsIGFuZCB0aGVuIHdlCmNhbiBnZXQgaXQgaW50byB0aGUgZm9ybWF0IHJlcXVpcmVkIGZvciBDS01SLgpgYGB7ciBrZWxwLWdlbm9zfQojIHdlIHdpbGwgdXNlIHRoaXMgc29tZSBtb3JlCmtnMiA8LSBub19oaV9taXNzZXJzICU+JSAKICBzZWxlY3QoTk1GU19ETkFfSUQsIGxvY3VzLCBhbGxlbGUpICU+JQogIG11dGF0ZShDaHJvbSA9ICJHVHNlcSIpICU+JSAKICBtdXRhdGUoUG9zID0gYXMuaW50ZWdlcihmYWN0b3IobG9jdXMsIGxldmVscyA9IHVuaXF1ZShsb2N1cykpKSkgJT4lCiAgcmVuYW1lKExvY3VzID0gbG9jdXMsCiAgICAgICAgIEFsbGVsZSA9IGFsbGVsZSkgJT4lCiAgc2VsZWN0KE5NRlNfRE5BX0lELCBDaHJvbSwgTG9jdXMsIFBvcywgQWxsZWxlKSAlPiUKICB1bmdyb3VwKCkKCiMgZ2V0IHRoZSBhbGxlbGUgZnJlcXMKa2dfY2ttcl9tYXJrZXJzIDwtIGtnMiAlPiUKICBmaWx0ZXIoIWlzLm5hKEFsbGVsZSkpICU+JSAjIGl0IGlzIHZpdGFsIHRvIGZpbHRlciBvdXQgdGhlIE5BcyBhdCB0aGlzIHN0YWdlCiAgZ3JvdXBfYnkoQ2hyb20sIExvY3VzLCBQb3MsIEFsbGVsZSkgJT4lCiAgc3VtbWFyaXNlKGNvdW50cyA9IG4oKSkgJT4lCiAgZ3JvdXBfYnkoTG9jdXMsIFBvcykgJT4lCiAgbXV0YXRlKEZyZXEgPSBjb3VudHMgLyBzdW0oY291bnRzKSkgJT4lCiAgc2VsZWN0KC1jb3VudHMpICU+JQogIG11dGF0ZShBbGxlSWR4ID0gMSwKICAgICAgICAgTG9jSWR4ID0gMSkgJT4lCiAgcmVpbmRleF9tYXJrZXJzKC4pCmBgYAoKCiMjIE1ha2luZyBnZW5vdHlwZSBtYXRyaWNlcyAmIGNoZWNraW5nIGZvciBkdXBsaWNhdGVkIHNhbXBsZXMKCmBgYHtyfQpub2Z1X2hhcHMgPC0gbm9faGlfbWlzc2VycyAlPiUKICBmaWx0ZXIoIWlzLm5hKGFsbGVsZSkpICU+JSAgIyBvbmNlIGFnYWluLCBpdCBpcyBjcml0aWNhbCB0byByZW1vdmUgdGhlc2UgYXQgdGhpcyBwb2ludAogIHNlbGVjdChOTUZTX0ROQV9JRCwgbG9jdXMsIGdlbmVfY29weSwgYWxsZWxlKSAlPiUKICByZW5hbWUoTG9jdXMgPSBsb2N1cywgQWxsZWxlID0gYWxsZWxlKQoKbm9mdV9pZHhfZnJhbWUgPC0ga2dfY2ttcl9tYXJrZXJzICU+JQogIHNlbGVjdChMb2N1cywgQWxsZWxlLCBMb2NJZHgsIEFsbGVJZHgpICU+JQogIGdyb3VwX2J5KExvY3VzKSAlPiUKICBtdXRhdGUoTnVtQSA9IG4oKSkgJT4lICAjIGdldCB0aGUgbnVtYmVyIG9mIGFsbGVsZXMgYXQgZWFjaCBsb2N1cwogIHVuZ3JvdXAoKSAlPiUKICBsZWZ0X2pvaW4obm9mdV9oYXBzLCAuKSAgJT4lICAjIGpvaW4gdGhlIGFsbGVfaWR4J3Mgb250byB0aGUgYWN0dWFsIGdlbm90eXBlIGRhdGEKICBzZWxlY3QoTk1GU19ETkFfSUQsIExvY3VzLCBnZW5lX2NvcHksIExvY0lkeCwgTnVtQSwgQWxsZUlkeCkgJT4lCiAgc3ByZWFkKGtleSA9IGdlbmVfY29weSwgdmFsdWUgPSBBbGxlSWR4KSAlPiUKICBtdXRhdGUoR2Vub0lkeCA9IGluZGV4X2FiKGEgPSBgMWAsIGIgPSBgMmAsIEEgPSBOdW1BKSkKCiMgbWFrZSBhIG1hdHJpeCBvZiBnZW5vdHlwZSBpbnRlZ2VycyAKd2lkZV9ub2Z1IDwtIG5vZnVfaWR4X2ZyYW1lICU+JQogIHNlbGVjdChOTUZTX0ROQV9JRCwgTG9jSWR4LCBHZW5vSWR4KSAlPiUKICBzcHJlYWQoZGF0YSA9IC4sIGtleSA9IExvY0lkeCwgdmFsdWUgPSBHZW5vSWR4KQpgYGAKCkRvbid0IGZvcmdldCB0byBzZXQgTkEncyB0byAwLCBhbmQgdGhlbiBkZWNyZWFzZSBlYWNoIHZhbHVlIGJ5IDE6CmBgYHtyIG1ha2UtbWF0fQpub2Z1X21hdCA8LSBhcy5tYXRyaXgod2lkZV9ub2Z1WywgLTFdKQpyb3duYW1lcyhub2Z1X21hdCkgPC0gd2lkZV9ub2Z1JE5NRlNfRE5BX0lECm5vZnVfbWF0W2lzLm5hKG5vZnVfbWF0KV0gPC0gMApub2Z1X21hdCA8LSBub2Z1X21hdCAtIDEKc3RvcmFnZS5tb2RlKG5vZnVfbWF0KSA8LSAgImludGVnZXIiCmBgYAoKIyMgTG9va2luZyBmb3IgZHVwbGljYXRlZCBzYW1wbGVzCgpMb29raW5nIGJhY2sgYXQgQW5keSBSYW1leSdzIFVTR1MgbWV0YSBkYXRhLCB0aGVyZSBhcmUgc29tZSBkdXBsaWNhdGUgc2FtcGxlcz8gT3IgYXQgbGVhc3QgYmlyZHMgdGhhdCBtYXkgaGF2ZSBiZWVuIHNhbXBsZWQgbXVsdGlwbGUgdGltZXMuIAoKV2UgY2FuIHF1aWNrbHkgbG9vayB0aHJvdWdoIHJvY2t5X21hdCBmb3IgcGFpcnMgb2YgaW5kaXZzIHdpdGggbG90cyBvZiBtYXRjaGluZyBnZW5vdHlwZXMuCmBgYHtyIGNoZWNrLWZvci1kdXBlc30KbWF0Y2hlcnMgPC0gcGFpcndpc2VfZ2Vub19pZChTID0gbm9mdV9tYXQsIG1heF9taXNzID0gMTIpICU+JQogIGFycmFuZ2UobnVtX21pc21hdGNoKSAlPiUKICBtdXRhdGUoTk1GU19ETkFfSURfMSA9IHJvd25hbWVzKG5vZnVfbWF0KVtpbmQxXSwKICAgICAgICAgTk1GU19ETkFfSURfMiA9IHJvd25hbWVzKG5vZnVfbWF0KVtpbmQyXSkKCm1hdGNoZXJzICU+JQogIGFycmFuZ2UoZGVzYyhudW1fbG9jKSkgIyBuZWVkIHRvIHVzZSB0aGUgZmlsdGVyZWQgZ2Vub3R5cGUgZGF0YSBmb3IgdGhpcy4uLgpgYGAKVGhlcmUgYXJlIDEzIHNhbXBsZXMgdGhhdCBtYXRjaCBwZXJmZWN0bHkuIEkgc2hvdWxkIGdldCByaWQgb2Ygb25lIG9mIGVhY2ggcGFpci4KCk5vdywgdG8gZGVhbCB3aXRoIHRoZXNlIGR1cGxpY2F0ZXMsIGhlcmUgaXMgd2hhdCBJIHdpbGwgZG86IHNpbmNlIGFsbCBvZiB0aGUgbWF0Y2hlcnMgaGF2ZSAwIG1pc21hdGNoZWQgbG9jaSwgd2UgYXJlIGdvaW5nIHRvIHRha2UganVzdCBvbmUgZnJvbSBlYWNoIHBhaXJzLiAKYGBge3IgdG9zcy1tYXRjaGVyc30KZGVmX3NhbWUgPC0gbWF0Y2hlcnMgJT4lCiAgZmlsdGVyKG51bV9taXNtYXRjaCA8PSAzKQoKdG9zc190aGVzZSA8LSBjKGludGVyc2VjdChkZWZfc2FtZSRpbmQxLCBkZWZfc2FtZSRpbmQyKSwgZGVmX3NhbWUkaW5kMikgJT4lCiAgdW5pcXVlKCkKCm5vZnVfbWF0X3Rvc3NlZCA8LSBub2Z1X21hdFstdG9zc190aGVzZSxdCmBgYAoKCmBgYHtyfQpub19oaV9taXNzZXJzX2NsZWFuIDwtIG5vX2hpX21pc3NlcnMgJT4lCiAgZmlsdGVyKE5NRlNfRE5BX0lEICVpbiUgcm93bmFtZXMobm9mdV9tYXRfdG9zc2VkKSkKCiMgc2F2ZSB0aGF0IGZvciB0aGUgYnljYXRjaCBhc3NpZ25tZW50cwp3cml0ZV9jc3Yobm9faGlfbWlzc2Vyc19jbGVhbiwgImNzdl9vdXRwdXRzL2Jhc2VsaW5lX25vX2hpX21pc3NlcnNfY2xlYW4uY3N2IikKYGBgCgoKCg==