vignettes/tutorial.Rmd
tutorial.Rmd
ADMIXTOOLS is a widely used software package for calculating admixture statistics and testing population admixture hypotheses.
A typical ADMIXTOOLS workflow generally involves a combination of sed
/awk
/shell scripting and manual editing to create text configuration files. These are then passed as command-line arguments to one of ADMIXTOOLS commands, and control how to run a particular analysis. The results are then redirected to another file, which has to be parsed by the user to extract values of interest, often using command-line utilities again or (worse) by manual copy-pasting. Finally, the processed results are analysed in R, Excel or another program.
This workflow can be a little cumbersome, especially if one wants to explore many hypotheses involving different combinations of populations. Most importantly, however, it makes it difficult to coduct reproducible research, as it is nearly impossible to construct fully automated “pipelines” that don’t require user intervention.
This R package makes it possible to perform all stages of ADMIXTOOLS analyses entirely from R, completely removing the need for “low level” configuration of individual ADMIXTOOLS programs.
Note that in order to use the admixr package, you need a working installation of ADMIXTOOLS! You can find installation instructions here. The software runs on Linux and macOS and these are the two systems that admixr is tested on.
Furthermore, you also need to make sure that R can find ADMIXTOOLS binaries on the $PATH
. You can achieve this by specifying PATH=<path to the location of ADMIXTOOLS programs>
in the .Renviron
file.
To install admixr you can simply run the following command in your R session:
install.packages("admixr")
Furthermore, if you want to follow the examples in this vignette, you will need the tidyverse collection of packages for data manipulation manipulation and plotting, which you can install with:
install.packages("tidyverse")
You definitely don’t need tidyverse for working with admixr but it really makes data manipulation and plotting things much easier. I recommend at least giving it a shot.
When everything is ready, you can run the following code to load both packages:
ADMIXTOOLS software uses a peculiar set of genetic file formats, which may seem strange if you are used to working with VCF files. However, the basic idea remains the same: we want to store and access SNP data (REF/ALT alleles) of a set of individuals at a defined set of genomic positions.
EIGENSTRAT datasets always contain three kinds of files:
ind
file - specifies a unique name, sex (optional - can be simply “U” for “undefined”) and label (such as population assignment) of each sample;snp
file - specifies the positions of SNPs, REF/ALT alleles etc.;geno
file - contains SNP data (one row per site, one character per sample) in a dense string-based format:
Therefore, a VCF file is essentially a combination of all three files in a single package.
Let’s first download a small testing SNP dataset using a built-in admixr function download_data()
. This function downloads the data into a temporary directory (you can specify the destination using its dirname
argument, in case you want to place it elsewhere).  Importantly, when specifying population/sample names in admixr functions, the information in the third column is what is used. For example, if you have individuals such as “French1”, “French2”, “French3” in the first column of an ind
file, all three sharing a “French” population label in the third column, specifying “French” in an admixr function will combine all three samples in a single population, instead of working with each individual separately.
The goal of admixr is to make ADMIXTOOLS analyses as trivial to run as possible, without having to worry about par/pop/left/right configuration files (as they are known in the jargon of ADMIXTOOLS) and other low-level details.
The only interface between you and ADMIXTOOLS is the following set of R functions:
Anything that would normally require dozens of lines of shell scripts can be often accomplished by running a single line of R code.
As we saw above, each EIGENSTRAT dataset has three components. The way this data is internally represented in admixr is using a small S3 R object created using the eigenstrat
constructor function. This function accepts the path and prefix of a trio of EIGENSTRAT snp/ind/geno files and returns an R object of the class EIGENSTRAT
:
snps <- eigenstrat(prefix)
snps #> EIGENSTRAT object #> ================= #> components: #> ind file: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/snps/snps.ind #> snp file: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/snps/snps.snp #> geno file: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/snps/snps.geno
This object encapsulates the paths to all three EIGENSTRAT components and makes it easy to pass the data to different admixr functions.
The following couple of sections describe how to use the admixr package using simple example analyses.
Let’s say we are interested in the following question: "Which populations today show evidence of Neanderthal admixture?
One way of looking at this is using the following D statistic: \[D(\textrm{present-day human W}, \textrm{African}, \textrm{Neanderthal}, \textrm{Chimp}).\]
\(D\) statistics are based on comparing the proportions of BABA and ABBA sites patterns observed in the data:
\[D = \frac{\textrm{# BABA sites - # ABBA sites}}{\textrm{# BABA sites + # ABBA sites}}.\]
Significant departure of \(D\) from zero indicates an excess of allele sharing between the first and the third population (positive \(D\)), or an excess of allele sharing between the second and the third population (negative \(D\)). If we get \(D\) that is not significantly different from 0, this suggests that the first and second populations form a clade, and don’t differ in the rate of allele sharing with the third population (this is the null hypothesis that the data is compared against).
Therefore, our \(D\) statistic above tests whether some modern humans today admixed with Neanderthals, which would increase their genetic affinity to this archaic group compared to Africans (whose ancestors never met Neanderthals).
Let’s save some population names first to make our code more concise:
pops <- c("French", "Sardinian", "Han", "Papuan", "Khomani_San", "Mbuti", "Dinka")
Using the admixr package we can then calculate our \(D\) statistic simply by running:
result <- d(W = pops, X = "Yoruba", Y = "Vindija", Z = "Chimp", data = snps)
The result is a following data frame:
head(result)
| W | X | Y | Z | D | stderr | Zscore | BABA | ABBA | nsnps |
|---|---|---|---|---|---|---|---|---|---|
| French | Yoruba | Vindija | Chimp | 0.0313 | 0.006933 | 4.510 | 15802 | 14844 | 487753 |
| Sardinian | Yoruba | Vindija | Chimp | 0.0287 | 0.006792 | 4.222 | 15729 | 14852 | 487646 |
| Han | Yoruba | Vindija | Chimp | 0.0278 | 0.006609 | 4.199 | 15780 | 14928 | 487925 |
| Papuan | Yoruba | Vindija | Chimp | 0.0457 | 0.006571 | 6.953 | 16131 | 14721 | 487694 |
| Khomani_San | Yoruba | Vindija | Chimp | 0.0066 | 0.006292 | 1.051 | 16168 | 15955 | 487564 |
| Mbuti | Yoruba | Vindija | Chimp | -0.0005 | 0.006345 | -0.074 | 15751 | 15766 | 487642 |
We can see that in addition to the specified population names, the output table contains additional columns:
D
- \(D\) statistic valuestderr
- standard error of the \(D\) statistic calculated using the block jackknifeZscore
- \(Z\)-zscore value (number of standard errors the \(D\) is from 0, i.e. how strongly do we reject the null hypothesis of no admixture)BABA
, ABBA
- counts of observed site patternsnsnps
- number of SNPs used for a given calculationWhile we could certainly make inferences by looking at the \(Z\)-scores, tables in general are not the best representation of this kind of data, especially as the number of samples increases. Instead, we can use the ggplot2
package to plot the results:
ggplot(result, aes(fct_reorder(W, D), D, color = abs(Zscore) > 2)) + geom_point() + geom_hline(yintercept = 0, linetype = 2) + geom_errorbar(aes(ymin = D - 2 * stderr, ymax = D + 2 * stderr))
(If you want to more know about data analysis using R, including plotting with ggplot2, I highly recommend this free book.)
We can see that the \(D\) values for Africans are not significantly different from 0, meaning that the data is consistent with the null hypothesis of no Neanderthal ancestry in Africans. On the other hand, the test rejects the null hypothesis for all non-Africans today, suggesting that Neanderthals admixed with the ancestors of present-day non-Africans.
An alternative way of addressing the previous question is to use the \(f_4\) statistic, which is very similar to \(D\) statistic and can be calculated as:
\[ f_4 = \frac{\textrm{# BABA sites - # ABBA sites}}{\textrm{# sites}}\]
Again, significant departure of \(f_4\) from 0 can be interpreted as evidence of gene flow.
To repeat the previous analysis using \(f_4\) statistic, we can run the function f4()
:
result <- f4(W = pops, X = "Yoruba", Y = "Vindija", Z = "Chimp", data = snps)
head(result)
| W | X | Y | Z | f4 | stderr | Zscore | BABA | ABBA | nsnps |
|---|---|---|---|---|---|---|---|---|---|
| French | Yoruba | Vindija | Chimp | 0.001965 | 0.000437 | 4.501 | 15802 | 14844 | 487753 |
| Sardinian | Yoruba | Vindija | Chimp | 0.001798 | 0.000427 | 4.209 | 15729 | 14852 | 487646 |
| Han | Yoruba | Vindija | Chimp | 0.001746 | 0.000418 | 4.178 | 15780 | 14928 | 487925 |
| Papuan | Yoruba | Vindija | Chimp | 0.002890 | 0.000417 | 6.924 | 16131 | 14721 | 487694 |
| Khomani_San | Yoruba | Vindija | Chimp | 0.000436 | 0.000415 | 1.051 | 16168 | 15955 | 487564 |
| Mbuti | Yoruba | Vindija | Chimp | -0.000030 | 0.000410 | -0.074 | 15751 | 15766 | 487642 |
By comparing this result to the \(D\) statistic analysis above, we can make the same conclusions.
You might be wondering why we have both \(f_4\) and \(D\) if they are so similar. The truth is that \(f_4\) is, among other things, directly informative about the amount of shared genetic drift (“branch length”) between pairs of populations, which is a very useful theoretical property. Other than that, it’s often a matter of personal preference and so admixr provides functions for calculating both.
Now we know that non-Africans today carry some Neanderthal ancestry. But what if we want to know how much Neanderthal ancestry they have? What proportion of their genomes is of Neanderthal origin? To answer questions like this, we can use the \(f_4\)-ratio statistic, which can be formulated in the following way (using a notation of Patterson et al., 2012, who formally described its properties).
\[f_4\textrm{-ratio} = \frac{f_4(A, O; X, C)}{f_4(A, O; B, C)}.\]
Using amidxr
, we can calculate \(f_4\)-ratios using the following code (X
being a vector of samples which we want to estimate the Neanderthal ancestry in):
result <- f4ratio(X = pops, A = "Altai", B = "Vindija", C = "Yoruba", O = "Chimp", data = snps)
The ancestry proportion (a number between 0 and 1) is given in the alpha
column:
head(result)
| A | B | X | C | O | alpha | stderr | Zscore |
|---|---|---|---|---|---|---|---|
| Altai | Vindija | French | Yoruba | Chimp | 0.023774 | 0.006176 | 3.850 |
| Altai | Vindija | Sardinian | Yoruba | Chimp | 0.024468 | 0.006071 | 4.031 |
| Altai | Vindija | Han | Yoruba | Chimp | 0.022117 | 0.005892 | 3.754 |
| Altai | Vindija | Papuan | Yoruba | Chimp | 0.037311 | 0.005812 | 6.420 |
| Altai | Vindija | Khomani_San | Yoruba | Chimp | 0.003909 | 0.005913 | 0.661 |
| Altai | Vindija | Mbuti | Yoruba | Chimp | 0.000319 | 0.005717 | 0.056 |
ggplot(result, aes(fct_reorder(X, alpha), alpha, color = abs(Zscore) > 2)) + geom_point() + geom_errorbar(aes(ymin = alpha - 2 * stderr, ymax = alpha + 2 * stderr)) + geom_hline(yintercept = 0, linetype = 2) + labs(y = "Neandertal ancestry proportion", x = "present-day individual")
We can make several observations:
The \(f_3\) statistic, also known as the 3-population statistic, is useful whenever we want to:
As an example, imagine we are interested in relative divergence times between pairs of present-day human populations, and want to know in which approximate order they split of from each other. To address this problem, we could use \(f_3\) statistic by fixing the \(C\) outgroup as San, and calculating pairwise \(f_3\) statistics between all present-day modern humans.
pops <- c("French", "Sardinian", "Han", "Papuan", "Mbuti", "Dinka", "Yoruba") result <- f3(A = pops, B = pops, C = "Khomani_San", data = snps)
head(result)
| A | B | C | f3 | stderr | Zscore | nsnps |
|---|---|---|---|---|---|---|
| French | French | Khomani_San | 0.000000 | -1.000000 | 0.000 | -1 |
| French | Sardinian | Khomani_San | 0.353447 | 0.012527 | 28.215 | 249760 |
| French | Han | Khomani_San | 0.316964 | 0.011914 | 26.604 | 253158 |
| French | Papuan | Khomani_San | 0.306962 | 0.011708 | 26.218 | 251648 |
| French | Mbuti | Khomani_San | 0.119283 | 0.008448 | 14.119 | 271501 |
| French | Dinka | Khomani_San | 0.190141 | 0.010049 | 18.922 | 276964 |
# sort the population labels according to an increasing f3 value relative to French ordered <- filter(result, A == "Mbuti", B != "Mbuti") %>% arrange(f3) %>% .[["B"]] %>% c("Mbuti") # plot heatmap of pairwise f3 values result %>% filter(A != B) %>% mutate(A = factor(A, levels = ordered), B = factor(B, levels = ordered)) %>% ggplot(aes(A, B)) + geom_tile(aes(fill = f3))
We can see that when we order the heatmap labels based on values of pairwise \(f_3\) statistics, the (already known) order of population splits pops up nicely (i.e. San separated first, followed by Mbuti, etc.).
Both qpWave and qpAdm can be though of as more complex and powerful extensions of the basic ideas behind a simple \(f_4\) statistic. Building upon the \(f_4\) theory and generalizing it, qpWave makes it possible to find the lowest number of “streams of ancestry” between two groups of populations that is consistent with the data. Extending the concept of \(f_4\) statistics even further, qpAdm allows to find the proportions of ancestry from a set of ancestral populations that contributed ancestry to our population of interest.
Unfortunately, both methods represent a rather advanced topic that still lacks proper documentation and beginner-friendly tutorials, and explaining them in detail is beyond the scope of this vignette. If you want to use them, it’s crucial that you read the official documentation decribing the basic ideas of both methods (distributed with ADMIXTOOLS), and that you read the relevant supplementary sections of papers published by David Reich’s group. At the very least, I recommend reading:
Note S6 of “Reconstructing Native American population history” by Reich et al. This paper first introduced the theoretical background of what later became qpWave.
Supplementary Information 10 of “Massive migration from the steppe was a source for Indo-European languages in Europe” by Haak et al., which gives a more consise overview of the qpWave method than S6 of Reich et al. 2012, and also introduces the qpAdm methodology for estimating admixture proportions.
A phenomenal description of qpAdm methodology and best practices by Harney et al. 2020 (see also PDF with practical guidelines. Both documents should really be a mandatory reading before doing any qpAdm analysis.
In the remainder of this section, I will assume that you are familiar with both methods, and will only explain how to use admixr for running them from R.
To run qpWave, you must provide a list of left and right populations (using the terminology of Haak et al. 2015 above). The aim of the method is to get an idea about the number of migration waves from right to left (with no back-migration from left to right!). This is done by estimating the rank of a matrix of all possible \(f_4\) statistics
\[f_4(\textrm{left}_1, \textrm{left}_i; \textrm{right}_1, \textrm{right}_i),\]
where \(\textrm{left}_1\) and \(\textrm{right}_1\) are some fixed populations and the \(i\) and \(j\) indices run over all other possible choices of populations.
As an example, let’s try to find the number of admixture waves from right = {Yoruba, Mbuti, Alta} into left = {French, Sardinian, Han} populations. We can do this using the function qpWave()
, setting its arguments appropriately:
result <- qpWave( left = c("French", "Sardinian", "Han"), right = c("Altai", "Yoruba", "Mbuti"), data = snps )
The qpWave()
function returns a data frame which shows the results of a series of matrix rank tests. The rank
column is the matrix rank tested, df
, chisq
and tail
give the degrees of freedom, \(\chi^2\) value and \(p\)-value for the comparison with the saturated model (the \(p\)-value then indicates which matrix rank is consistent with the data - see example below), and dfdiff
, chisqdiff
and taildiff
give the same, but always comparing a model to the model with one rank less.
result
| rank | df | chisq | tail | dfdiff | chisqdiff | taildiff |
|---|---|---|---|---|---|---|
| 0 | 4 | 1.756 | 0.7805106 | 0 | 0.000 | 1.0000000 |
| 1 | 1 | 0.192 | 0.6614968 | 3 | 1.564 | 0.6674998 |
| 2 | 0 | 0.000 | 1.0000000 | 1 | 0.192 | 0.6614968 |
In this example, we see that matrix \(r = 0\) cannot be rejected (tail
\(p\)-value = 0.78). Because Reich et al. 2012 showed that \(r + 1 \le n\), where \(n\) is the number of admixture waves, we can interpret this as left populations having at least \(n = 1\) streams of ancestry from the set of right populations. In this case, the most likely explanation is Neandertal admixture into non-Africans today.
Now, what happens if we add Papuans to the left group?
result <- qpWave( left = c("Papuan", "French", "Sardinian", "Han"), right = c("Altai", "Yoruba", "Mbuti"), data = snps )
result
| rank | df | chisq | tail | dfdiff | chisqdiff | taildiff |
|---|---|---|---|---|---|---|
| 0 | 6 | 29.143 | 0.0000572 | 0 | 0.000 | 1.0000000 |
| 1 | 2 | 0.601 | 0.7403643 | 4 | 28.542 | 0.0000097 |
| 2 | 0 | 0.000 | 1.0000000 | 2 | 0.601 | 0.7403643 |
We can now clearly reject rank \(r = 0\), but we see that the data is consistent with rank \(r = 1\), meaning that there must have been at least \(n = 2\) streams of ancestry from right to left populations (\(r + 1 \le n\)). Because this happened after we introduced Papuans to the left set, this could indicate a separate pulse of archaic introgression into Papuans, which is not surprising given what we know about significantly more archaic ancestry in Papuans than in any other present-day population.
The qpAdm method can be used to find, for a given target population, the proportions of ancestry coming from a set of source populations. Importantly, since we often lack accurate representatives of the true ancestral populations, we can use a set of reference populations instead, under a crucial assumption that the references set is phylogenetically closer to true source populations than to a set of specified outgroups. For example, coming back to our example of estimating the proportions of Neandertal ancestry in people today, we could define:
Having defined all three population sets, we can run qpAdm with:
result <- qpAdm( target = c("Sardinian", "Han", "French"), sources = c("Vindija", "Yoruba"), outgroups = c("Chimp", "Denisova", "Altai"), data = snps )
The qpAdm()
function has an argument details
(default TRUE) which makes the function return a list of three elements:
proportions
- data frame with admixture proportions - this is what we mostly care about;ranks
- results of rank tests performed by qpWave - these evaluate how well does the assumed traget-sources-outgroups population model match the data;subsets
- results of the “all subsets” analysis (see the documentation for more details.If details
is set to FALSE
, only the proportions
components is returned by the qpAdm()
function.
Let’s start with the ranks
element:
result$ranks
| target | rank | df | chisq | tail | dfdiff | chisqdiff | taildiff |
|---|---|---|---|---|---|---|---|
| Sardinian | 1 | 1 | 0.005 | 0.9447256 | 3 | -0.005 | 1.0000000 |
| Sardinian | 2 | 0 | 0.000 | 1.0000000 | 1 | 0.005 | 0.9447256 |
| Han | 1 | 1 | 2.096 | 0.1476639 | 3 | -2.096 | 1.0000000 |
| Han | 2 | 0 | 0.000 | 1.0000000 | 1 | 2.096 | 0.1476639 |
| French | 1 | 1 | 3.745 | 0.0529780 | 3 | -3.745 | 1.0000000 |
| French | 2 | 0 | 0.000 | 1.0000000 | 1 | 3.745 | 0.0529780 |
The row with rank = 1 represents a qpWave test with all \(n\) source populations set as the left set and all outgroups as the right set. This test evaluates whether the ancestral populations are descended from \(n\) independent streams of ancestry. In our case, \(n = 2\) (Mbuti and Vindija), which means that the data would have to be consistent with rank \(r = 1\) to satisfy the inequality \(r + 1 \le n\) proved by Reich et al., 2012. We see that this is true for all three target populations (\(p\)-value > 0.05 for all targets), and the simple model of Neandertal admixture thus seems to be reasonably consistent with the data.
The rank = 2 row represents a qpWave test after adding a target population to the left group together with the sources. This test makes sure that including the target population does not increase the rank of the \(f_4\) matrix, meaning that the target can be really modelled as a mixture of ancestries from the sources. If the \(p\)-values turn out to be very low, this indicates that the assumed model does not fit the data and that a part of the ancestry in a target possibly cannot be traced to any of the sources. In our case, however, all rank = 2 test \(p\)-values are not significant, and we can be reasonably sure that the target samples can be fully modelled as a mixtures of all specified references.
The most important element of a qpAdm output is in the $proportions
component. This contains admixture proportion estimates from all specified sources, the p-values of each model (remember, low p-values/significance means the model is rejected!) as well as standard errors for those proportions using a block jackknife:
result$proportions
| target | Vindija | Yoruba | stderr_Vindija | stderr_Yoruba | nsnps_used | nsnps_target | pvalue |
|---|---|---|---|---|---|---|---|
| Sardinian | 0.024 | 0.976 | 0.006 | 0.006 | 5e+05 | 499314 | 0.944726 |
| Han | 0.020 | 0.980 | 0.006 | 0.006 | 5e+05 | 499654 | 0.147664 |
| French | 0.022 | 0.978 | 0.006 | 0.006 | 5e+05 | 499434 | 0.052978 |
If we compare this result to the \(f_4\)-ratio values calculated above, we see that the qpAdm estimates are very close to what we got earlier.
The third element in the list of results shows the outcome of an “all subsets” analysis, which involves testing all subsets of potential source populations. Each 1 in the “pattern” column means that the proportion of ancestry from that particular source population (in the order as specified by the user) was forced to 0.0.
result$subsets
| target | pattern | wt | dof | chisq | tail | Vindija | Yoruba |
|---|---|---|---|---|---|---|---|
| Sardinian | 00 | 0 | 1 | 0.005 | 0.9447260 | 0.024 | 0.976 |
| Sardinian | 01 | 1 | 2 | 15795.353 | 0.0000000 | 1.000 | 0.000 |
| Sardinian | 10 | 1 | 2 | 16.295 | 0.0002895 | 0.000 | 1.000 |
| Han | 00 | 0 | 1 | 2.096 | 0.1476640 | 0.020 | 0.980 |
| Han | 01 | 1 | 2 | 14809.047 | 0.0000000 | 1.000 | 0.000 |
| Han | 10 | 1 | 2 | 14.239 | 0.0008092 | 0.000 | 1.000 |
| French | 00 | 0 | 1 | 3.745 | 0.0529780 | 0.022 | 0.978 |
| French | 01 | 1 | 2 | 15326.662 | 0.0000000 | 1.000 | 0.000 |
| French | 10 | 1 | 2 | 16.102 | 0.0003188 | 0.000 | 1.000 |
New feature: There is a new function called qpAdm_rotation()
which allows exhaustive exploration of many qpAdm models in paralle. For more information please see the vignette _“Fitting qpAdm models with a ‘rotation’ strategy_”.
What we’ve been doing so far was calculating statistics for individual samples. However, it is often useful to treat multiple samples as a single group or population. admixr provides a function called relabel()
that does just that.
Here is an example: let’s say we want to run a similar analysis to the one described in the \(D\) statistic section, but we want to treat Europeans, Africans and archaics as combined populations, and not as separate individuals. But the ind
file that we have does not contain grouped labels - each sample stands on its own:
Chimp U Chimp
Mbuti U Mbuti
Yoruba U Yoruba
Khomani_San U Khomani_San
Han U Han
Dinka U Dinka
Sardinian U Sardinian
Papuan U Papuan
French U French
Vindija U Vindija
Altai U Altai
Denisova U Denisova
To merge several individual samples under a combined label we can call relabel()
like this:
# paths to the original ind file and a new modified ind file, which will # contain merged population labels modif_snps <- relabel( snps, European = c("French", "Sardinian"), African = c("Dinka", "Yoruba", "Mbuti", "Khomani_San"), Archaic = c("Vindija", "Altai", "Denisova") )
modif_snps #> EIGENSTRAT object #> ================= #> components: #> ind file: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/snps/snps.ind #> snp file: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/snps/snps.snp #> geno file: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/snps/snps.geno #> #> modifiers: #> groups: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/file79e927f7c386.ind
We can see that the function relabel
returned a modified EIGENSTRAT
object, which contains a new item in the “modifiers” section - the path to a new ind file.  For example:
result <- d(W = "European", X = "African", Y = "Archaic", Z = "Chimp", data = modif_snps)
Here is the result, showing again Europeans show genetic affinity to archaic humans compared to Africans today:
head(result)
| W | X | Y | Z | D | stderr | Zscore | BABA | ABBA | nsnps |
|---|---|---|---|---|---|---|---|---|---|
| European | African | Archaic | Chimp | 0.0225 | 0.004404 | 5.117 | 15487 | 14805 | 489003 |
Note that the d()
function correctly picks up the “group modifier” ind
file from the provided EIGENSTRAT
object and uses it in place of the original ind
file.
The count_snps
function can be useful for quality control, weighting of admixture statistics (\(D\), \(f_4\), etc.) in regression analyses etc. There are two optional arguments:
prop
- changes whether to report SNP counts or proportions (set to FALSE
by default),missing
- controls whether to count missing SNPs instead of present SNPs (set to FALSE
by default).For each sample, count the SNPs present in that sample:
count_snps(snps)
| id | sex | label | present |
|---|---|---|---|
| Chimp | U | Chimp | 491273 |
| Mbuti | U | Mbuti | 499334 |
| Yoruba | U | Yoruba | 499246 |
| Khomani_San | U | Khomani_San | 499250 |
| Han | U | Han | 499654 |
| Dinka | U | Dinka | 499362 |
| Sardinian | U | Sardinian | 499314 |
| Papuan | U | Papuan | 499377 |
| French | U | French | 499434 |
| Vindija | U | Vindija | 497544 |
| Altai | U | Altai | 497729 |
| Denisova | U | Denisova | 497398 |
It is quite common to repeat a particular analysis only on a subset of the genome (such as intergenic sites, etc). However, EIGENSTRAT is a rather obscure file format which is generally not supported by standard bioinformatics tools. Luckily, admixr includes a function filter_bed()
that takes an EIGENSTRAT
object and a BED file as its inputs and produces a new object that contains a modifier called “excluded”, linking to a snp file with coordinates of sites that did not pass the filtering and will be excluded from later analyses.
# BED file contains regions to keep in an analysis new_snps <- filter_bed(snps, bed) # BED file contains regions to remove from an analysis new_snps <- filter_bed(snps, bed, remove = TRUE)
new_snps #> EIGENSTRAT object #> ================= #> components: #> ind file: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/snps/snps.ind #> snp file: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/snps/snps.snp #> geno file: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/snps/snps.geno #> #> modifiers: #> excluded sites: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/file79e96c229701.snp #> (SNPs excluded: 100000, SNPs remaining: 400000)
If we want to run the whole analysis in a single pipeline, we can use the %>%
pipe operator and do the following:
(The %>%
operator takes what is on its left side and puts it as a first argument of a function on the right side. While it takes some time to get used to, it is very useful in longer multi-step “pipelines” because it makes more pipelines much more readable. In fact, the resulting code often reads almost like English! The %>%
pipe is automatically imported when you load the tidyverse
library, and you can read about it more here.)
snps %>% filter_bed("regions.bed") %>% d(W = "French", X = "Mbuti", Y = "Vindija", Z = "Chimp")
This is because in the formal definitions of admixr function, data =
is always the argument, so we don’t have to specify it manually.
Important: The filter_bed()
function makes it very easy to do filtering without worrying about locations of intermediate files, but it is important to keep in mind that the function still creates temporary files under the hood. If you plan to run many independent calculations on a filtered subset of the data, it’s better to save the new EIGENSTRAT
object to a variable first and re-use the same object multiple times, rather than running the whole pipeline for each analysis separately (which would create new copies of intermediate files for each iteration).
In the field of ancient DNA, we often need to repeat an analysis on a subset of data that is less likely to be influenced by ancient DNA damage, to verify that our results are not caused by artifacts in the data (due to biochemical properties of DNA degradation, ancient DNA damage will lead to an increase in C→T and G→A substitutions). Using a similar method described in the BED filtering section above, we can use the transversions_only()
function to generate a snp file with positions that carry transitions (C→T and G→A sites):
new_snps <- transversions_only(snps) # perform the calculation only on transversions d(W = "French", X = "Dinka", Y = "Altai", Z = "Chimp", data = new_snps)
Again, we could combine several filtering steps into one pipeline:
snps %>% # take the original data filter_bed("regions.bed", remove = TRUE) %>% # remove sites not in specified regions transversions_only() %>% # remove potential false SNPs due to aDNA damage d(W = "French", X = "Dinka", Y = "Altai", Z = "Chimp") # calculate D on the filtered dataset
Another useful data processing function is merge_eigenstrat()
. This function takes two EIGENSTRAT datasets and merges them, producing a union of samples and intersection of SNPs from both of them and returning a new EIGENSTRAT
object.
# this is just an example code - it will not run unless you specify the paths merged <- merge_eigenstrat( merged = <"prefix of the merged dataset"> a = first_EIGENSTRAT_object, b = second_EIGENSTRAT_object )
The goal of admixr is to abstract away all the low-level technical workings of ADMIXTOOLS. As we saw in the examples above, it achieves this by doing all the dirty work of parsing the output files generated by ADMIXTOOLS, presenting the user with convenient R data structures.
Nevertheless, admixr cannot (yet) parse all the information generated by ADMIXTOOLS utilities. For calculating D statistics, \(f_4\) statistics and doing other simple analyses, a single data frame is usually all the user needs. However, sometimes things go wrong - it turns out data was in a wrong format, or ADMIXTOOLS crashes with an error or the results simply look suspicious. Furthermore, some complex commands (such as qpAdm) are not yet completely implemented in admixr. In any of these cases, it is still useful (and important!) to examine the log outputs.
Each admixr result, such as the data frame object below, contains a (hidden!) attribute which carries the complete log file associate with that particular analysis:
dres <- d(W = c("French", "Han", "Dinka"), X = "Yoruba", Y = "Vindija", Z = "Chimp", data = snps) dres #> # A tibble: 3 x 10 #> W X Y Z D stderr Zscore BABA ABBA nsnps #> <chr> <chr> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> #> 1 French Yoruba Vindija Chimp 0.0313 0.00693 4.51 15802 14844 487753 #> 2 Han Yoruba Vindija Chimp 0.0278 0.00661 4.20 15780 14928 487925 #> 3 Dinka Yoruba Vindija Chimp -0.0009 0.00612 -0.151 15131 15159 487667
We can examine the full log of this D statistic run with the loginfo
function. Without any further arguments, this function writes the complete log on the screen. Note that we run this function on the data frame object although it might appear that it does not contain any other information (it does but it’s normally hidden as most of the time it’s not useful):
loginfo(dres) #> Full output log of qpDstat: #> =================================================== #> #> qpDstat: parameter file: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/qpDstat__2016673616.par #> ### THE INPUT PARAMETERS #> ##PARAMETER NAME: VALUE #> genotypename: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/snps/snps.geno #> snpname: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/snps/snps.snp #> indivname: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/snps/snps.ind #> popfilename: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/qpDstat__2016673616.pop #> printsd: YES #> ## qpDstat version: 970 #> number of quadruples 3 #> 0 French 1 #> 1 Han 1 #> 2 Dinka 1 #> 3 Yoruba 1 #> 4 Vindija 1 #> 5 Chimp 1 #> jackknife block size: 0.050 #> snps: 500000 indivs: 6 #> number of blocks for block jackknife: 547 #> nrows, ncols: 6 500000 #> result: French Yoruba Vindija Chimp 0.0313 0.006933 4.510 15802 14844 487753 #> result: Han Yoruba Vindija Chimp 0.0278 0.006609 4.199 15780 14928 487925 #> result: Dinka Yoruba Vindija Chimp -0.0009 0.006124 -0.151 15131 15159 487667 #> ##end of qpDstat: 3.807 seconds cpu 0.000 Mbytes in use
The situation is a bit more complex for qpAdm which can be evaluated for multiple target populations/samples which are then analyzed separately, with each analysis having its own log file:
qpadm_res <- qpAdm( target = c("Sardinian", "Han"), sources = c("Vindija", "Yoruba"), outgroups = c("Chimp", "Denisova", "Altai"), data = snps ) qpadm_res #> $proportions #> # A tibble: 2 x 8 #> target Vindija Yoruba stderr_Vindija stderr_Yoruba nsnps_used nsnps_target #> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> #> 1 Sardi… 0.024 0.976 0.006 0.006<PHONE_NUMBER>14 #> 2 Han 0.02 0.98 0.006 0.006<PHONE_NUMBER>54 #> # … with 1 more variable: pvalue <dbl> #> #> $ranks #> # A tibble: 4 x 8 #> target rank df chisq tail dfdiff chisqdiff taildiff #> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> #> 1 Sardinian 1 1 0.005 0.945 3 -0.005 1 #> 2 Sardinian 2 0 0 1 1 0.005 0.945 #> 3 Han 1 1 2.10 0.148 3 -2.10 1 #> 4 Han 2 0 0 1 1 2.10 0.148 #> #> $subsets #> # A tibble: 6 x 8 #> target pattern wt dof chisq tail Vindija Yoruba #> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> #> 1 Sardinian 00 0 1 0.005 0.945 0.024 0.976 #> 2 Sardinian 01 1 2 15795. 0 1 0 #> 3 Sardinian 10 1 2 16.3 0.000289 0 1 #> 4 Han 00 0 1 2.10 0.148 0.02 0.98 #> 5 Han 01 1 2 14809. 0 1 0 #> 6 Han 10 1 2 14.2 0.000809 0 1
In this case, we could either print all log information (one log output for each target we specified in the qpAdm call) by runningloginfo(qpadm_res)
or, perhaps more useful, specify which target’s log file we want to examine:
loginfo(qpadm_res, target = "Han") #> Full output log of qpAdm for target 'Han': #> =================================================== #> #> qpAdm: parameter file: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/qpAdm__1149292693.par #> ### THE INPUT PARAMETERS #> ##PARAMETER NAME: VALUE #> genotypename: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/snps/snps.geno #> snpname: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/snps/snps.snp #> indivname: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/snps/snps.ind #> popleft: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/qpAdm__1149292693.popleft #> popright: /var/folders/t7/9gjtb6m92flbnp930618vt3r0000gn/T//RtmpveArl0/qpAdm__1149292693.popright #> allsnps: YES #> summary: YES #> details: YES #> ## qpAdm version: 1201 #> seed: 436781117 #> #> left pops: #> Han #> Vindija #> Yoruba #> #> right pops: #> Chimp #> Denisova #> Altai #> #> 0 Han 1 #> 1 Vindija 1 #> 2 Yoruba 1 #> 3 Chimp 1 #> 4 Denisova 1 #> 5 Altai 1 #> jackknife block size: 0.050 #> snps: 500000 indivs: 6 #> number of blocks for block jackknife: 547 #> ## ncols: 500000 #> coverage: Han 499654 #> coverage: Vindija 497544 #> coverage: Yoruba 499246 #> coverage: Chimp 491273 #> coverage: Denisova 497398 #> coverage: Altai 497729 #> Effective number of blocks: 489.823 #> numsnps used: 500000 #> codimension 1 #> f4info: #> f4rank: 1 dof: 1 chisq: 2.096 tail: 0.147663931 dofdiff: 3 chisqdiff: -2.096 taildiff: 1 #> B: #> scale 1.000 #> Denisova 0.475 #> Altai 1.332 #> A: #> scale 27.655 #> Vindija 1.414 #> Yoruba -0.029 #> #> #> full rank #> f4info: #> f4rank: 2 dof: 0 chisq: 0.000 tail: 1 dofdiff: 1 chisqdiff: 2.096 taildiff: 0.147663931 #> B: #> scale 19.608 759.097 #> Denisova 0.470 -0.810 #> Altai 1.334 -1.160 #> A: #> scale 1.414 1.414 #> Vindija 1.414 0.000 #> Yoruba 0.000 1.414 #> #> #> best coefficients: 0.020 0.980 #> Jackknife mean: 0.020445464 0.979554536 #> std. errors: 0.006 0.006 #> #> error covariance (* 1,000,000) #> 33 -33 #> -33 33 #> #> #> summ: Han 2 0.147664 0.020 0.980 33 -33 33 #> #> fixed pat wt dof chisq tail prob #> 00 0 1 2.096 0.147664 0.020 0.980 #> 01 1 2 14809.047 0 1.000 0.000 #> 10 1 2 14.239 0.000809235 0.000 1.000 #> best pat: 00 0.147664 - - #> best pat: 10 0.000809235 chi(nested): 12.143 p-value for nested model: 0.000492824 #> #> coeffs: 0.020 0.980 #> #> ## dscore:: f_4(Base, Fit, Rbase, right2) #> ## genstat:: f_4(Base, Fit, right1, right2) #> #> details: Vindija Denisova 0.023947 40.866731 #> details: Yoruba Denisova -0.001067 -2.572745 #> dscore: Denisova f4: -0.000556 Z: -1.350781 #> #> details: Vindija Altai 0.068033 120.675791 #> details: Yoruba Altai -0.001528 -3.697347 #> dscore: Altai f4: -0.000107 Z: -0.260309 #> #> gendstat: Chimp Denisova -1.351 #> gendstat: Chimp Altai -0.260 #> gendstat: Denisova Altai 1.118 #> #> ##end of qpAdm: 3.633 seconds cpu 0.000 Mbytes in use
Finally, we might want to keep some of the log files around for future reference, further debugging or sharing with others. For these purposes, we can run the loginfo
function with the argument save = TRUE
, potentially also specifying the output directory (dir = "."
by default) or a prefix of the output file(s) (default prefix is simply the name of the admixr command which generated the log file).
For example, the following will save the result of our qpAdm analysis of Neandertal ancestry to a file with a prefix “qpAdm_Neandertal_ancestry”, but just the one for the Sardinian individual:
loginfo(qpadm_res, target = "Sardinian", save = TRUE, prefix = "qpAdm_Neandertal_ancestry")
This will create a file in the current directory named qpAdm_Neandertal_ancestry_Sardinian.txt
. Note that when you call loginfo()
on a qpAdm run, the name of the target is always added to the end of the output log file.
This work is licensed under a Creative Commons Attribution 4.0 International License.


--- layout: page title: Lean Prover Zulip Chat Archive permalink: archive/113488general/07618preordercomap.html ---
Do we already have the function which given a map X -> Y and a preorder on Y, produces a preorder on X? I searched for preorder.comap and found nothing. Is that what it's supposed to be called?
You'd think this just gets synthesised automatically when preorder
is defined. :grinning:
(More seriously --- how close could we get to that happening? Can we construct derive handlers so a new structure could just be labelled with @[derive comap]
and that definition is added?)
I think it's called preorder.lift


California school shooting suspect, 16, dies in hospital a day after gun rampage
Nathaniel Berhow is believed to have opened fire at Saugus High School in Santa Clarita, California, using a handgun hidden in his backpack - before turning the weapon om himself
by Joshua Taylor, https://www.facebook.com/dailymirror/A 16-year-old boy suspected of killing two classmates in a school shooting has died in hospital.
Nathaniel Berhow is said to have shot himself in the head after killing the students, aged 16 and 14, at Saugus High School in Santa Clarita, California.
The suspect died this afternoon with his mother next to him after his life support machine was switched off, according to police.
The shooter pulled a handgun from his backpack and emptied it in a matter of seconds as the school day began.
The gun rampage on Thursday - Berhow's 16th birthday - left two students dead and three more injured.
Captain Kent Wegener, of the Los Angeles County Sheriff's Department, said the entire incident, captured on videotape, lasted just 16 seconds as the young man stood in one spot and fired at one student after another.
"From right where he was standing, he doesn't chase anybody, he fires from where he is until he shoots himself," Wegener said.
The scene awas reminiscent of other mass shootings at US schools, including Marjory Stoneman Douglas High School in Parkland, Florida, where a former student with an assault gun killed 17 people in February 2018.
It was the 85th incident of gunfire at a school this year, according to Everytown, a gun control advocacy group.
Authorities have no clues about what sparked this latest bloodshed.
"We did not find any manifesto, any diary that spelled it out," Los Angeles County Sheriff Alex Villanueva said at a briefing.
Two girls, aged 14 and 15, remain in hospital in stable conditions. A boy, 14, who was also injured, has been released.
The FBI is assisting in the investigation and police do not know yet where the suspect got the weapon.
Katie Holt, a fine arts teacher, pushed the piano in the music room against the classroom door to create a barricade.
"She turned off the lights... She locked all the doors. She pushed the piano to the door, totally barricading," said Eddie Mendoza, 17 in a phone interview, his voice still breaking from emotion and grief.
One student was heard saying "I think I got shot" as a group locked themselves inside a choir classroom.


Magister Ladd:
Some heroes reject destiny, and some have it thrust upon them... For instance, Raul, the Unification King, was chosen by the Goddess of the Lake, and there was Emperor Shunaiman, who received power from Einhasad... Don't you agree that they were chosen by the gods?
Those who are chosen have no use for Mimir's Elixir.
But the Star of Destiny... that's another story...
Do you want to make Mimir's Elixir? How about it? Can you handle the power?
"Yes."


joe starred jsmith/wordpad 1 hours ago
ansuz forked derp/crypto to ansuz/crypto 2 hours ago
pro starred twitter/bootstrap 3 hours ago
ninja starred jgthms/bulma 3 hours ago
username starred name/space 19 hours ago
username starred name/space 19 hours ago
username starred name/space 19 hours ago
username starred name/space 19 hours ago
username starred name/space 19 hours ago
username starred name/space 19 hours ago
username starred name/space 19 hours ago
username starred name/space 19 hours ago
username starred name/space 19 hours ago
username starred name/space 19 hours ago
username starred name/space 19 hours ago
username starred name/space 19 hours ago


This is a strong element. This element portrays information that the writer wants the reader to catch their attention.
This is a bolded element. This element depicts important information in text.
This is an encoded element. This element emphasises a word.
This is an italicized element. This element italicizes a word.
This is a quotation
element. This element puts quotation marks around a word.
This is an source code enclosed
element. This element puts code into the html page.
This is an keyboard encolsed element.
This is an enclosed above subscript element. This element places smaler text above the rest of the text to work as a note to refer the reader to some other part of the text.
This is an enclosed below subscript element. This element places smaller text below the rest of the text to work as another way to refer the reader to another part of the text.
This is an incline logical span divisiont element. This element is used to group a few elements together in an organized format.
I did not see much change from normal text in the span element. It just seemed to look like normal text with no formating unlike the other inline elements.


The {@link oaj.annotation.Swap#template() @Swap(template)} annotation allows you to associate
arbitrary contextual strings with swaps.
The primary purpose is for providing template names, such as for Apache FreeMarker, therefore
the name 'template'.
However, the usage of the string is open-ended.
For example, you could pair a template string like so:
| @Swap(impl=FreeMarkerSwap.class, template="MyPojo.div.ftl")
| public class MyPojo {}
The implementation of the FreeMarker swap would look something like this:
| // Our templated swap class.
| public class FreeMarkerSwap extends ObjectSwap<Object,Reader> {
|
| public MediaType[] forMediaTypes() {
| // Make sure this only applies to the HTML serializer.
| return MediaType.forStrings("*/html");
| }
|
| public Reader swap(BeanSession session, Object object, String template) throws Exception {
| // Call some method that uses FreeMarker to convert 'o' to raw HTML using
| // the 'MyPojo.div.ftl' template.
| return getFreeMarkerReader(template, object);
| }
| }


Class HcsDid
- java.lang.Object
-
- com.hedera.hashgraph.identity.hcs.did.HcsDid
-
-
Field Summary
Fields Modifier and Type Field Description static DidSyntax.Method
DID_METHOD
-
Constructor Summary
Constructors Constructor Description HcsDid(java.lang.String network, com.hedera.hashgraph.sdk.PrivateKey privateDidRootKey, com.hedera.hashgraph.sdk.FileId addressBookFileId, com.hedera.hashgraph.sdk.TopicId didTopicId)
Creates a DID instance with private DID root key.HcsDid(java.lang.String network, com.hedera.hashgraph.sdk.PublicKey didRootKey, com.hedera.hashgraph.sdk.FileId addressBookFileId)
Creates a DID instance without topic ID specification.HcsDid(java.lang.String network, com.hedera.hashgraph.sdk.PublicKey didRootKey, com.hedera.hashgraph.sdk.FileId addressBookFileId, com.hedera.hashgraph.sdk.TopicId didTopicId)
Creates a DID instance.HcsDid(java.lang.String network, java.lang.String idString, com.hedera.hashgraph.sdk.FileId addressBookFileId, com.hedera.hashgraph.sdk.TopicId didTopicId)
Creates a DID instance.
-
Method Summary
Modifier and Type Method Description static HcsDid
fromString(java.lang.String didString)
Converts a Hedera DID string intoHcsDid
object.DidDocumentBase
generateDidDocument()
Generates a basic DID document for this identity.DidDocumentBase
generateDidDocument(com.hedera.hashgraph.sdk.PublicKey didRootKey)
Generates DID document base from the given DID and its root key.static com.hedera.hashgraph.sdk.PrivateKey
generateDidRootKey()
Generates a random DID root key.com.hedera.hashgraph.sdk.FileId
getAddressBookFileId()
com.hedera.hashgraph.sdk.TopicId
getDidTopicId()
java.lang.String
getIdString()
DidSyntax.Method
getMethod()
Returns a DID method of this DID.java.lang.String
getNetwork()
Returns aString
for which the DID is created.java.util.Optional<com.hedera.hashgraph.sdk.PrivateKey>
getPrivateDidRootKey()
Returns a private key of DID root key.static java.lang.String
publicKeyToIdString(com.hedera.hashgraph.sdk.PublicKey didRootKey)
Constructs an id-string of a DID from a given public key.java.lang.String
toDid()
Converts theHederaDid
object into a DID string.java.lang.String
toString()
-
-
-
Field Detail
-
DID_METHOD
public static final DidSyntax.Method DID_METHOD
-
-
Constructor Detail
-
HcsDid
public HcsDid(java.lang.String network, com.hedera.hashgraph.sdk.PublicKey didRootKey, com.hedera.hashgraph.sdk.FileId addressBookFileId, com.hedera.hashgraph.sdk.TopicId didTopicId)
Creates a DID instance.- Parameters:
network
- The Hedera DID network.didRootKey
- The public key from which DID is derived.addressBookFileId
- The appent's address bookFileId
didTopicId
- The appnet's DID topic ID.
-
HcsDid
public HcsDid(java.lang.String network, com.hedera.hashgraph.sdk.PrivateKey privateDidRootKey, com.hedera.hashgraph.sdk.FileId addressBookFileId, com.hedera.hashgraph.sdk.TopicId didTopicId)
Creates a DID instance with private DID root key.- Parameters:
network
- The Hedera DID network.privateDidRootKey
- The private DID root key.addressBookFileId
- The appent's address bookFileId
didTopicId
- The appnet's DID topic ID.
-
HcsDid
public HcsDid(java.lang.String network, com.hedera.hashgraph.sdk.PublicKey didRootKey, com.hedera.hashgraph.sdk.FileId addressBookFileId)
Creates a DID instance without topic ID specification.- Parameters:
network
- The Hedera DID network.didRootKey
- The public key from which DID is derived.addressBookFileId
- The appent's address bookFileId
-
HcsDid
public HcsDid(java.lang.String network, java.lang.String idString, com.hedera.hashgraph.sdk.FileId addressBookFileId, com.hedera.hashgraph.sdk.TopicId didTopicId)
Creates a DID instance.- Parameters:
network
- The Hedera DID network.idString
- The id-string of a DID.addressBookFileId
- The appent's address bookFileId
didTopicId
- The appnet's DID topic ID.
-
-
Method Detail
-
fromString
public static HcsDid fromString(java.lang.String didString)
Converts a Hedera DID string intoHcsDid
object.- Parameters:
didString
- A Hedera DID string.- Returns:
HcsDid
object derived from the given Hedera DID string.
-
generateDidRootKey
public static com.hedera.hashgraph.sdk.PrivateKey generateDidRootKey()
Generates a random DID root key.- Returns:
- A private key of generated public DID root key.
-
publicKeyToIdString
public static java.lang.String publicKeyToIdString(com.hedera.hashgraph.sdk.PublicKey didRootKey)
Constructs an id-string of a DID from a given public key.- Parameters:
didRootKey
- Public Key from which the DID is created.- Returns:
- The id-string of a DID that is a Base58-encoded SHA-256 hash of a given public key.
-
toDid
public java.lang.String toDid()
Description copied from interface:HederaDid
Converts theHederaDid
object into a DID string.
-
generateDidDocument
public DidDocumentBase generateDidDocument()
Description copied from interface:HederaDid
Generates a basic DID document for this identity.- Specified by:
generateDidDocument
in interfaceHederaDid
- Returns:
- The DID document as a
DidDocumentBase
.
-
generateDidDocument
public DidDocumentBase generateDidDocument(com.hedera.hashgraph.sdk.PublicKey didRootKey)
Generates DID document base from the given DID and its root key.- Parameters:
didRootKey
- Public key used to build this DID.- Returns:
- The DID document base.
- Throws:
java.lang.IllegalArgumentException
- In case given DID root key does not match this DID.
-
getNetwork
public java.lang.String getNetwork()
Description copied from interface:HederaDid
Returns aString
for which the DID is created.- Specified by:
getNetwork
in interfaceHederaDid
- Returns:
- A Hedera network for which the DID is created.
-
getMethod
public DidSyntax.Method getMethod()
Description copied from interface:HederaDid
Returns a DID method of this DID.
-
toString
public java.lang.String toString()
- Overrides:
toString
in classjava.lang.Object
-
getDidTopicId
public com.hedera.hashgraph.sdk.TopicId getDidTopicId()
-
getAddressBookFileId
public com.hedera.hashgraph.sdk.FileId getAddressBookFileId()
-
getIdString
public java.lang.String getIdString()
-
getPrivateDidRootKey
public java.util.Optional<com.hedera.hashgraph.sdk.PrivateKey> getPrivateDidRootKey()
Returns a private key of DID root key. This is only available if it was provided duringHcsDid
construction.- Returns:
- The private key of DID root key.
-
-


185 Stevens Way AC 101
Paul G. Allen Center Box 352350
Seattle, WA 98185-2350
<EMAIL_ADDRESS>Publications
Kyle Rector, Lauren R. Milne, Richard E. Ladner, Batya Friedman, Julie A. Kientz. 2015. Exploring the Opportunities and Challenges with Exercise Technologies for People who are Bind or Low-Vision. In Proc. ASSETS ’15. ACM, New York, NY, USA.
Catherine M. Baker, Lauren R. Milne, and Richard E. Ladner. 2015. StructJumper: A Tool to Help Blind Programmers Navigate and Understand the Structure of Code. In Proceedings of the ACM Conference on human factors in computing systems (CHI ‘15). ACM, New York, NY, USA.
Lauren R. Milne, Cynthia L. Bennett, Shiri Azenkot, and Richard E. Ladner. 2014. BraillePlay: Educational Smartphone Games for Blind Children. In Proceedings of the 16th International ACM SIGACCESS Conference on Computers and Accessibility (ASSETS '14). ACM, New York, NY, USA.
Best Student Paper:
Catherine M. Baker, Lauren R. Milne, Jeffrey Scofield, Cynthia L. Bennett, and Richard E. Ladner. 2014. Tactile Graphics with a Voice: Using QR Codes to Access Text in Tactile Graphics. In Proceedings of the 16th International ACM SIGACCESS Conference on Computers and Accessibility (ASSETS '14). ACM, New York, NY, USA.
Lauren R. Milne, Cynthia L. Bennett, and Richard E. Ladner. 2014. The Accessibility of Mobile Health Sensors for Blind Users. In 29th Annual International & Persons with Disabilities Conference (CSUN ’14).
Conference Posters
Catherine M. Baker, Lauren R. Milne, Jeffrey Scofield, Cynthia L. Bennett, and Richard E. Ladner. 2014. Tactile Graphics with a Voice Demonstration. In Proceedings of the 16th International ACM SIGACCESS Conference on Computers and Accessibility (ASSETS '14). ACM, New York, NY, USA.
Lauren R. Milne, Cynthia L. Bennett, and Richard E. Ladner. 2013. VBGhost: a Braille-Based Educational Smartphone Game for Children. In Proceedings of the 15th International ACM SIGACCESS Conference on Computers and Accessibility (ASSETS '13). ACM, New York, NY, USA, Article 75, 2 pages.


Chapter 14
THE WARNERS AND RELIGION
The unknown Warners who crossed
the
A ritual drinking feast with alcohol provided magical knowledge. This magical knowledge is typically associated with the quest for good fortune. The religious ideas of the Anglo-Saxons made allowance for the worship of many gods. The names of their many Gods became a part of our calendar. Their language, Angl-ish, became our English language.
Sunne, the sun god became our Sunday; Mona, the moon god Monday; Twi, the god of warfare and battle Tuesday; Woden - ancestor of people, leader of the hunt and carries off the dead Wednesday; Thunor, the god of thunder, who rules the storms and sky Thursday; and Frige, the goddess of love and the wife of Woden became our Friday.
In 597, Pope Gregory sent
In 780 Charlemagne decreed the death penalty for all Saxons
who failed to be baptized, who failed to keep Christian festivals, and who
cremated their dead. By the sword and the cross, Charlemagne became master of
As the Anglo-Saxon culture developed, those first known Warners were probably Catholic. The power of king and church were combined in mutual support of one another. The national throne and the church attempted to control what the people were to believe and how they were to worship. The growing awareness of individual rights and natural curiosity led to opinions and discoveries in nature that challenged both church teachings and tradition. Early attempts to reform the church led to the Reformation and several new Protestant denominations.
The Puritans were Calvinist in their theology and believed in the harsh judgment of God. They embraced predestination and rejected free will. Quakers on the other hand saw God as a Spirit who lives in the heart of every human being and offers direct and personal guidance. Personal experience with God is taken to be more authoritative than the Bible. The Warners were Quakers in the 1500s until Captain William Warner (7) chose to join with the Puritans.
The Puritans came to
Captain William Warner and family traveled with a group of
Quakers when they came to
The Puritans and Quakers differed on their belief about the
relationship between man and God. The Puritans believed human beings were
depraved and unable to commune with God. Only by grace (that is, by God's
decision) could any individual achieve salvation. There was nothing man could
do to achieve it.
Quakers believed that Christ was like a bridge between God and man. They
believed the Sermon on the Mount was to be taken literally and, by obeying
Jesus' commands, that individuals could earn salvation. What the Quakers and
the Puritans had in common was the belief that, in the end, the vast majority
of people would not be saved.
Members of the Quaker movement became trusted and respected for their honesty and were described as living the faith that others distorted. They sought to live in peace with everyone, including the Native Americans, and practiced religious toleration within their communities.
We can imagine the Captain s grandson, William Warner (9), sitting in the usual Quaker Meeting where someone had just spoken about the treatment of slaves. When they gathered for their Meetings, they sat in the circle. All people were considered as created equal and could speak but not interrupt another. There would be meditation and consideration of what the Lord would want. During a silence or meditation, a revelation was thought as direct from God. Feeling such was enough to make one s body shake hence the name Quaker. I wonder how William reacted when he felt that God wanted him to forget the wealth of his investment and release his nine slaves.
Six generations of Warners were
Quakers. The son of Mordecai (12) and Jane, John Lewis (13), was the last. He had grown up a
Quaker but did not join a Meeting after leaving
Many of the children of John and Charity became active in
the
The church near the house MH Warner built in
When RJ wrote about religion and the relatives back in
The older boys had grown up in a green world a world of clean air, creeks, orchards, woods, meadows, pastures and dogs, sheep and cattle, horses and farmyards. They read selections from the great historians, orators, novelists, essayists, and poets. The McGuffy Reader was a literary storehouse for family reading and a portable library for ambitious youths in a nearly bookless country.
|
Sixth McGuffey Reader used by both MH Warner and RJ Warner Now in possession of Elvira Warner Covey Published 1879 |
RJ wrote about religion:
A heaven to gain and a hell to shun, was the way the pronouncement of the
evangelists who introduced my generation of our family to religion, at the
beginning of the twentieth century. The broad road ahead led to the destruction
and eternal damnation. The narrow road led to a life of eternal bliss. Each road
had a guide to assist us in our passage through this life. The devil would direct
us down the broad road and would introduce us to all the sins that humanity
enjoyed, while the Son of God , would lead us down the narrow road, through
many trials and temptations, but to a victorious future. We were urged to choose
a leader that would make us a winner.
The message of Christ s saving death was offered as a choice for empowerment to live life for eternal good. The transforming and positive power of this belief appears regularly in the conversion stories of many individuals who find joy in God s love and grace. RJ had a framed scripture that always hung above the headboard of his bed, Lo I am with you always. Matthew 28:20.
A few years after he was married and bought the house in New Rome, RJ took his family to a church just down the street. During the service he noticed that every time a certain man got up to speak they immediately sang a hymn. After a while an usher came over to him and said they were having a problem and apologized. It seems the man wanted to say something about the theory of evolution and the congregation did not want to hear him. The next Sunday RJ took his family to another church and that may be why we were raised as Methodists.
|
While I was growing up we always attended Sunday school and church. Sometimes we would get up early, go part way to visit a relative and attend a a church along the way. Mom was a stickler for perfect attendance (fifty Sundays a year) and made sure we earned our Sunday school attendance pins. My brother, R. David had pins for thirteen years. Mine, shown here is for ten years. Dad always taught a Sunday school class. His choice of Bible was the King James with the Thompson chain reference. Mom had several Bibles and usually added the latest translation to the bookshelf beside her chair. She started the junior church when the teaching method was telling Bible stories using a flannel board. The Bible characters came in a book, were colored with crayons, cut out with scissors and placed on the board as the story was told. |
Dad gave me a Bible when I left home. He wrote this inscription inside the front cover.
Dear Bob, Christmas 1943
This is the best book I have ever known. It has helped me and it will help you wherever
you go. Study its pages and you will find secret strength for all the trials of life: guidance
for all the dark pathways that might lead to wrong, and you will develop character that
will stand every test man knows. When the years of life ripen into old age, if you have
accepted Christ, you can look back and see that you have lived a good and useful life.
Dad
While I was growing up I never felt that I was pressured to get religion but rather I felt that my parents had discovered something that I should catch up with. Like my Dad I attended revivals at other denominations and came to my own conclusions. The Christian story has much to offer: power to overcome wrong choices, release from the burden of guilt and reproach, service to one s neighbor, and sensitivity to the suffering of the least. It offers a sense of confidence that life has meaning and purpose.
The MH Warner family was a religious family. Perhaps having parents in different denominations developed a greater awareness of what was really important in matters of faith. They had attained a confidence that goes with grasping truths that are above the meaning of words. In something like the Quaker tradition they were confident of their ability to ascertain what was of God and what was not. They left the hills unencumbered with many of the restrictions imposed by some believers. Many members of the following generations are active in various churches and several are preachers.
Go to Chapter:
1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11,
12, 13,
14, 15,
16, 17


Planet Sick-Boy: http://www.sick-boy.com "We Put the SIN in Cinema"
© Copyright 2001 Planet Sick-Boy. All Rights Reserved.
As cheesy as the original Planet of the Apes film is, it sets a pretty high cinematic bar. It'd be easy to outdo the special effects and best the acting, but tackling the infamous Statue of Liberty ending is an entirely different story. Duplicating or removing it would be a cop-out, and changing it virtually guarantees the new version won't be any good (it'd be like making Rosebud a woman in a Citizen Kane remake).
While the new ending is a little clunky (and a lot predictable), director Tim Burton does, for the most part, a pretty decent job with his "revisitation" of Apes. It's probably the most un-Burton-esque film he has made, relying more on special effects, makeup and costume than his typically stunning visual presentation. Apes is a rare remake that is shorter than its predecessor, and one that doesn't cast actors from the original...because most of them are dead (aside from Charlton Heston, who is dead only between his ears, and makes a brief appearance here). There are enough differences to make this film unique (for better or for worse), but it's also similar enough to avoid sending diehard Ape fans into nerd comas.
Apes opens in 2029 inside a research station just outside Saturn. Mark Wahlberg (The Perfect Storm) plays Leo Davidson, a pilot chomping at the bit to fly into an electrical storm (didn't Wahlberg learn to stay away from storms last summer?), but a space chimp is tapped to make the journey before any human lives are risked. When the monkey's pod disappears, Leo takes off after it, despite the orders to stop from his commanding officer. Sure enough, Leo's pod vanishes into a Contact-like wormhole that spits him out several hundred years later as he crashes into a strange planet ruled by apes.
Like the first film, there are humans on the planet, but they can talk and have more interaction with their monkey overlords, who act a lot more simian than they did in the original. Some humans are slaves, some are kept as pets and some have even been molded into butlers, but the apes are firmly in control of things despite their apparent lack of ability to invent stuff beyond rudimentary weapons and armor. There are no guns, no doctors and generally no knowledge of what lies over the mountain in The Forbidden Zone.
While Wahlberg's Leo is clearly crafted in the same mold as Charlton Heston's Colonel Taylor, the only ape character that makes the evolutionary leap from 1968 to 2001 is Ari, played by Fight Club's Helena Bonham Carter (her ancestor was Dr. Zira, which is almost Ari spelled backward). Ari, the rebellious daughter of a senator and the ex-lover of the evil General Thade (Tim Roth, Lucky Numbers), leads the Human Rights Faction and, for some reason, develops an instant bond with Leo.
With a few exceptions, the second act follows the original film pretty closely and finds Leo, his hot human squeeze (Estella Warren, Driven), Ari and a handful of others trying to make it into The Forbidden Zone while being chased by the gorilla army (led by Michael Clarke Duncan, See Spot Run). But as things degenerate into Braveheart (Apeheart?) or, worse yet, The 13th Warrior, it becomes clear this film is working its way to a more upbeat resolution, which kind of blows, and, coming from Burton, is just plain disappointing. Yeah, there's a twist, but you can see it coming from a mile away (assuming you know how the first film ended).
One major improvement in the film is the appearance of the apes. Rick Baker's (Dr. Seuss' How the Grinch Stole Christmas) creatures look amazing and can show much more emotion in their faces.  Another minor effects complaint is the over-the-top jumping skills of the chimps. I can overlook the fact that it's preposterous, but the editing of these scenes leaves a little to be desired. You never once see a clean landing, instead getting a jumpy edit. It's like throwing a dummy off a cliff and then having an actor pop into the shot from the bottom.
Apes, based on Pierre Boulle's novel, only scratches the surface of the complex racial issues and battles over religion versus science that made the first film such a treat. The screenplay, from William Broyles, Jr. (Cast Away), Lawrence Konner and Mark D. Rosenthal (Mighty Joe Young), does riff on a few of the original's more popular lines, and includes more children and teenage characters, both ape and human. What's more amazing is the inclusion of a swimming scene for Warren, who got to breaststroke in Driven, too (she used to be a synchronized swimmer).
Maybe it's because half of the actors were buried under latex (or because there's no Johnny Depp), but Apes seems like Burton's least actor-driven film. There's no weird cross-species chemistry between Wahlberg and Bonham Carter, but there isn't any between him and Warren, either. Roth sounds nothing like himself but becomes a wonderful snarling antagonist. If chimps had moustaches, he'd be twirling his. Danny Elfman's (Proof of Life) score reaches its peak during the opening credits, never factors into the rest of the film, and never equals the unmistakable power of Jerry Goldsmith's original.
1:47 - PG-13 for some sequences of action/violence
========== X-RAMR-ID: 28991 X-Language: en X-RT-ReviewID: 237328 X-RT-TitleID: 1108704 X-RT-SourceID: 595 X-RT-AuthorID: 1146
The review above was posted to the
rec.arts.movies.reviews newsgroup (de.rec.film.kritiken for German reviews).
The Internet Movie Database accepts no responsibility for the contents of the
review and has no editorial control. Unless stated otherwise, the copyright
belongs to the author.
Please direct comments/criticisms of the review to relevant newsgroups.
Broken URLs inthe reviews are the responsibility of the author.
The formatting of the review is likely to differ from the original due
to ASCII to HTML conversion.
Related links: index of all rec.arts.movies.reviews reviews


Learning resources
Internet of Things
Internet of Things 6 course specialization (Coursera)
Software Architecture for the Internet of Things (Coursera)
Cybersecurity and the Internet of Things (Coursera)
A developer's guide to the Internet of Things (by IBM)
The Raspberry Pi Platform and Python Programming for the Raspberry Pi (Coursera)


Code of Conduct
v 2019.02
Our Pledge
In the interest of fostering an open and welcoming environment, we as community members pledge to making participation in our community a harassment-free experience for everyone, regardless of age, body size, disability, ethnicity, sex characteristics, gender identity and expression, level of experience, education, socio-economic status, nationality,
This Code of Conduct applies to all of the community’s activities, including the Slack team, its code repositories, and its face-to-face meetups.
Our Standards
Examples of behaviour that contributes to creating a positive environment include:
- Using welcoming and inclusive language.
- Being respectful of differing viewpoints and experiences.
- Gracefully accepting constructive criticism.
- Focusing on what is best for the community.
- Showing empathy towards other community members.
Examples of unacceptable behaviour by participants include:
- The use of sexualized language or imagery and unwelcome sexual attention or advances.
- Trolling, insulting/derogatory comments, and personal or political attacks.
- Public or private harassment.
- Publishing others’ private information, such as a physical or electronic address, without explicit permission.
- Other conduct which could reasonably be considered inappropriate
Our Responsibilities
Community maintainers have the right and responsibility to remove, edit, or reject comments, code, edits, issues, and other contributions that are not aligned to this Code of Conduct, or to ban temporarily or permanently any contributor for other behaviours that they deem inappropriate, threatening, offensive, or harmful.
Scope
This Code of Conduct applies both within community spaces and in public space when an individual is representing the project or its community. Examples of representing a project or community include using an official project e-mail address, posting via an official social media account, or acting as an appointed representative at an online or offline event. Representation of the community may be further defined and clarified by community maintainers.
Enforcement
Instances of abusive, harassing, or otherwise unacceptable behaviour may be reported by contacting the maintainers at<EMAIL_ADDRESS>All complaints will be reviewed and investigated and will result in a response that is deemed necessary and appropriate to the circumstances. The community maintainers are obligated to maintain confidentiality with regard to the reporter of an incident. Further details of specific enforcement policies may be posted separately.
Community maintainers who do not follow or enforce the Code of Conduct in good faith may face temporary or permanent repercussions as determined by other members of the community.
This Code of Conduct is adapted from the Contributor Covenant, version 1.4, available here. For answers to common questions about this code of conduct, see here.


Deploying {{ project["name"] }}
Do not close this page while deploying.
Deploy
Waiting for deploy.
Back


mhmm [options] <motif file>
This program creates motif-based hidden Markov models (HMMs) of families of related biosequences. The program takes as input a set of DNA or protein motif models and produces as output a single HMM containing the given motifs. The program can produce three types of models: linear models, in which the motifs are arranged like beads on a string, completely connected models, which allow for repetitions of motifs and for motifs to appear in any order, and star models. Mhmm writes its output in a format readable by the other MEME Suite programs, mhmms and mhmmscan.
Three types of models may be produced. A linear motif-based HMM consists of a sequence of motif models, each separated by one or more tied insert states that represent the spacer region between motifs. A completely connected model, on the other hand, includes transitions from the end of each motif to the beginning of every other motif in the model (with a spacer model along each transition). This more general topology allows for motifs that are repeated, deleted or shuffled. A star model is also available. By default, the program produces a linear model.
Transition probabilities among motifs are derived from the motif occurrence information in the given MEME file. For the completely connected topology, this information is derived from all of the motif occurrences. For the linear topology, the information is derived only from the best-scoring sequence. Alternatively, the order and spacing of motifs within a linear model may be specified via the --order option.
A file containing MEME formatted motifs. Outputs from MEME and DREME are supported along with minimal MEME format for which there are conversion scripts available to support other formats. Input motifs that are likely to appear in the sequences.
Writes a motif-based hidden Markov model in MHMM format to standard output.
| Option | Parameter | Description | Default Behavior | ||||||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Input/Output | |||||||||||||||
| --motif | motif # | This option (which may be repeated) allows the user to select a specific motif for inclusion in the HMM. The specified motif number corresponds to the motif index in the MEME file. | |||||||||||||
| --nmotifs | n | This option is similar to the -motif option, except that it tells mhmm to use the first n motifs in the given MEME file. | All motifs are included. | ||||||||||||
| --ethresh | ev | This option sets an E-value threshold for inclusion of motifs in the model. Motifs with E-values < ev are ignored. | All motifs are included. | ||||||||||||
| --lowcomp | threshold | Eliminate low-complexity motifs from the model. Motif
complexity is the average K-L distance between the "motif
background distribution" and each column of the motif. The
motif background is just the average distribution of all the
columns. The K-L distance, which measures the difference
between two distributions, p and f, is the same as the
information content:
p1.log(p1/f1) +
p2.log(p2/f2) + ... +
pn.log(pn/fn)
This value increases with increasing complexity.
|
|||||||||||||
| --order | string | This option instructs mhmm to build an HMM with linear
topology. The given string specifies the order and spacing of
the motifs within the model, and has the format
"s=m=s=m=…=s=m=s ", where "s" is the length of a spacer
between motifs, and "m" is a motif ID. Thus, for example, the
string "34=3=17=2=5 " specifies a two-motif linear model, with
motifs 3 and 2 separated by 17 letters and flanked by 34 letters
and 5 letters on the left and right. If the MEME file contains
motif occurrences on both strands, then the motif IDs in the
order string should be preceded by "+" or "-" indicating the
strandedness of the motif. If this option is specified then the
--type option must either be
linear or not set. |
|||||||||||||
| --type | linear|complete|star |
|
A linear HMM is created. | ||||||||||||
| --pthresh | p-value | This option sets a p-value threshold for inclusion of motif occurrences in the transition probability matrix used to construct the HMM. | |||||||||||||
| --nspacer | num | Normally the mhmm program models each spacer using a single insert state. The distribution of spacer lengths produced by a single insert stage is exponential in form. A more reasonable distribution would be a bell-shaped curve such as a Gaussian. Modeling the length distribution explicitly is computationally expensive; however, as a Gaussian distribution can only be approximated using multiple insert states to represent a single spacer region. The --nspacer option specifies the number of insert states used to represent each spacer. Note that this Gaussian approximation is only effective in conjunction with total probability training and scoring. | |||||||||||||
| --transpseudo | pseudocount | Specify the value of the pseudocount used in converting transition counts to transition probabilities. | |||||||||||||
| --spacerpseudo | pseudocount | Specify the value of the pseudocount used in converting transition counts to spacer self-loop probabilities. | |||||||||||||
| --description | text | Specify descriptive text to be stored as a comment in the MHMM file. | |||||||||||||
| --fim | A free-insertion module (FIM) is an insert state with 1.0 probability of self-transition and 1.0 probability of exit transition. Thus, traversing such a state has zero transition cost. Specifying this option causes all spacers to be represented using FIMs. | ||||||||||||||
| --keep-unused | By default, mhmm remove from the transition probability matrix all inter-motif transitions that are not observed in the data. This option allows those transitions to be retained. This option is only relevant if the model has a completely connected topology. | ||||||||||||||
| --noheader | Do not put a header on the output file. | ||||||||||||||
| --noparams | Do not list the parameters at the end of the output. | ||||||||||||||
| --notime | Do not print the running time and host name at the end of the output. | ||||||||||||||
| --quiet | Combine the previous three flags and set verbosity to 1. |


View data anywhere around the world in minutes. It is secure and transparent at your ease.
Accessible, trusted and reliable source of monitoring and tracking CRYPTOCURRENCIES!!


Reactor Stream is a suite of Reactive Extensions operators and more implementing Reactive Streams. More detailed documentation is available on the reference guide.


Chapter 5 Average and SD
5.1 The Intuition
If we took a histogram and tried to describe it to someone else without showing it to them, the most important pieces of information are usually the average and standard deviation or SD.
We might describe the variable this way: “The values are about __________, give or take ________ or so.” We can think of the first blank as the average and the second blank as the SD.
- The average describes where the histogram is positioned along the left-right axis.
- The SD describes the width (or “spread” or “dispersion”) of the histogram.
Note that the average and SD are not the only way to describe the “location” and “width” (“dispersion,” “spread”, or “scale”) of the histogram, but they are the most common.
Inspect the collection of 20 numbers below. Fill in the blanks: “The values are about __________, give or take ________ or so.”
| 63 | 72 | 78 | 92 | 81 |
| 62 | 75 | 74 | 81 | 88 |
| 68 | 79 | 97 | 88 | 82 |
| 79 | 72 | 98 | 85 | 86 |
If we had to describe these data, we might say that our variable is “about 80, give or take 10 or so.” We can also visually these two summaries with a histogram.
While this variable has a particular average (about 80), we can imagine shifting it left or right. The figure below shows some possible shifts. We could shift it way to the left, so that it’s “about 40” or a little bit to the right so that it’s “about 90.”
We can also imagine increasing the SD (more spread) or decreasing the SD (less spread). The figure below shows some possible changes in the SD. In each case, the “give or take” number is changing.
Exercise 5.1 For the following nine lists, choose whether the average is 0, 5, or 10 and the SD is 2, 10, or 20. You don’t need to do any calculation.
- -22, -10, 35, 7, -29, 14, 9, -3, 23, 27
- -14, 17, 5, 8, 16, 12, 14, 16, 3, 24
- 10, -13, 6, 6, 4, -1, -8, 19, 6, 21
- 27, 37, -11, -23, 24, 15, -10, 2, 36, 2
- 5, 6, 6, 2, 3, 3, 9, 4, 6, 6
- -12, -3, 9, 11, -1, 2, -7, 3, 16, -18
- 12, 12, 9, 10, 12, 12, 6, 10, 11, 7
- 0, -2, 0, 3, 3, 0, -2, -3, -1, 2
- -3, 28, 12, 7, 25, 3, -38, -25, 7, -16
Solution
| List | Mean | SD |
|---|---|---|
| 1 | 5 | 20 |
| 2 | 10 | 10 |
| 3 | 5 | 10 |
| 4 | 10 | 20 |
| 5 | 5 | 2 |
| 6 | 0 | 10 |
| 7 | 10 | 2 |
| 8 | 0 | 2 |
| 9 | 0 | 20 |
Exercise 5.2 For the following nine histograms, choose whether the average is 0, 5, or 10 and the SD is 2, 10, or 20. You don’t need to do any calculations.
Solution
| Histogram | Mean | SD |
|---|---|---|
| Histogram 1 | 9.97 | 20.18 |
| Histogram 2 | 10.07 | 1.91 |
| Histogram 3 | 9.94 | 9.97 |
| Histogram 4 | 5.05 | 2.10 |
| Histogram 5 | 5.07 | 10.10 |
| Histogram 6 | 0.08 | 20.04 |
| Histogram 7 | 4.94 | 20.09 |
| Histogram 8 | -0.02 | 1.97 |
| Histogram 9 | 0.06 | 9.97 |
5.2 The Measures
5.2.1 The Average
The most common measure of the location of a variable is the average.1 Suppose we have a variable (a list of numbers) \(X = \{x_1, x_2, ..., x_n\}\).
\[\begin{equation} \text{average} = \dfrac{\text{the sum of the list}}{\text{the number of entries in the list}} = \dfrac{\sum_{i = 1}^n x_i}{n} \nonumber \end{equation}\]
The average is easy to compute and easy to work with mathematically.2
Unfortunately, the average doesn’t have an easy interpretation. The best interpretation, in my mind, is as the balance-point for the data. If we imagine the left-right axis as a teeter-totter and stack the data along the beam according to their values, then the average is the position of the fulcrum that would balance the data-filled beam.
Exercise 5.3 Compute the average of the following lists of numbers. Round to two decimals placed as necessary.
- {0, 1}
- {0, 0, 1}
- {0, 1, 1}
- {0, 10}
- {0, 1, 2}
- {3, 6, 7, 1}
- {8, -4, 16, 4, 7}
- {7, -4, 7, 7, -4}
Exercise 5.5 To find the percent of observations that fall into a particular category, I count the number of observations in that category and divide that by the total number. That number lies between 0 and one. We refer to that number as a “proportion.” To convert the proportion to a percent, we just multiply the proportion times 100%.
Suppose I asked 10 survey respondents whether they voted in the previous presidential election and obtained the following responses: {Yes, No, No, No, Yes, Yes, No, No, No, No}. What percent said “Yes”? What if I changed the “Yes” responses to the number 1 and the “No” responses to the number 0, then found the average? Is that related to the percent that said “Yes”? Is counting the “Yes” responses like changing them to ones, changing everything else to 0, and adding up the 0s and 1s?Solution
If the average of the ten numbers is 50, then the sum of the ten numbers must be 500. If the sum of the ten numbers is 500, then the sum of the eleven numbers is 600. Since the sum of the eleven numbers is 600, the average is 600/11 = 54.55.
5.2.2 The Median
When we’re describing the location of a variable (or, visually, of the histogram), we typically use the average. However, some researchers use the median. I introduce the median as an alternative to the average for two reasons. First, it helps us better understand the average conceptually. Second, you might encounter it, so it’s good to understand how the median works.
The median is the (or one of the) numbers that divides the variable (or, visually, the histogram) into two equal pieces.
For example, the number 3 divides the variable \(X = \{1, 2, 2, 3, 4, 6, 8\}\) into two equal pieces. The three values 1, 2, and 2 fall below the number 3. The three values 4, 6, and 8 fall above 3. Because 3 splits the variable into two equal pieces, it is the median.
If you order the variable from highest to lowest, you can think of the median as the “middle” value.
Notice, though, that things can run amok when we have an even number of observations in our variable. Consider the variable \(X = \{1, 2, 2, 4, 6, 8\}\). There is no “middle” value here.
Notice that the values 3, 2.7, and 3.76421 all split the variable into equal pieces (1, 2, and 2 fall below each of these values; 4, 6, and 8 fall above). This means that there are mmany values that split the variable into equal pieces. But it isn’t very helpful—we need a single median to summarize the variable.
When the number of observations in the variable is even, we just find the two middle values and split the difference. For the variable \(X = \{1, 2, 2, 4, 6, 8\}\), the numbers 2 and 4 are the two middle numbers. To “split the difference,” we just average them. The average of 2 and 4 is (2 + 4)/2 = 3. In this case, then, the median is three.
In summary, we find the median with the following algorithm:
- Order the values from lowest to highest.
- If the number of values is odd, then just find the middle value.
- If the number of values is even, then find the two middle values and split the difference (or average them).
Exercise 5.7 Find the average and median of the following variables. (I already ordered them from highest to lowest.)
- \(X = \{1, 3, 4, 6, 8\}\)
- \(X = \{1, 3, 4, 6, 80\}\)
- \(X = \{1, 3, 4, 6, 800\}\)
- \(X = \{1, 3, 4, 600, 800\}\)
Do you notice something interesting? Explain.
What is the average and median for \(X = \{1, 3, 400, 600, 800\}\)? What changed?The previous exercise illustrates something important about the average and the median. When a variable has a histogram with a longer tail in one direction, that long tail pulls the average away from the median in the direction of the tail. This happens because a adding few unusually large values can change average quite a lot, while changing the median very little.
Remember that the average is the balance point for the histogram. As a tail stretches out in one direction, those data exert greater leverage and can pull the average in that direction.
When a histogram has a longer tail to the right, we say that histogram is “right-skewed” or “skewed to the right.” When a histogram has a longer tail to the left, we say that histogram is “left-skewed” or “skewed to the left.”
The figure below shows examples of a histogram with no skew, one that is skewed to the right, and one that is skewed to the left. It shows the average and median for each. Notice how the long tail pulls the average away from the median toward the long tail.
## Rows: 6
## Columns: 4
## $ group <chr> "No Skew", "No Skew", "Skewed to the Left", "Skewed to the Left…
## $ name <chr> "Average", "Median", "Average", "Median", "Average", "Median"
## $ value <dbl> 51.08652, 51.05596, 75.87061, 78.36764, 19.71639, 17.28709
## $ y <dbl> 0, 0, 0, 0, 0, 0
Exercise 5.8 For each histogram below, roughly estimate the average. Next, decide whether the histogram has no skew, is skewed to the right, or is skewed to the left. Then decide whether the median is less than, greater than, or about equal to the average.
Solution
- no skew; median about the average
- no skew; median about the average
- skewed to the right; median below the average
- skewed to the left; median avove the average
- no skew; median about the average
5.2.3 The Standard Deviation
The most common measure of spread is the standard deviation (SD). The intuition is subtle, so let’s look a a simple example. Rember, our goal is a “give-or-take” number.
Suppose we have a list of numbers \(X = \{1, 2, 3, 4, 5\}\). The average of this list is 3, so we can compute the deviation from average or \(d\) for each value.
\[\begin{equation} \text{deviation from average} = d = \text{value} - \text{average} \nonumber \end{equation}\]
In this case, \(d = \{-2, -1, 0, 1, 2\}\).
We want to use these deviations to find a give-or-take number.
Here’s an initial idea. Just take the absolute values \(|d| = \{2, 1, 0, 1, 2\}\). These tell us how far each entry falls away from the average. Then we could average the absolute devations to find how far a typical entries falls away from the average of the list. In this case, we get 1.2. This is reasonable approach and we’ll refer to it as the “average absolute devation” or a.a.d. (It turns out that the a.a.d. isn’t a common quantity, so I don’t elevate it with an all-caps acronym.)
The a.a.d. has one big problem–it uses an absolute value. This introduces some computational and mathematical difficulties.
So let’s do something similar. Rather than take the absolute value, let’s square the deviations, take the average, and then undo the square at the end, so that \(\text{SD} = \sqrt{\text{avg}(d^2)}\).
Sometimes taking the (3) square root of (2) the average of (1) the squares is called the RMS. In this case, the RMS of the deviations from the average is the SD, so that
\[\begin{equation} \text{SD} = \sqrt{\text{avg}(d^2)} = \sqrt{\dfrac{(x_i - \text{avg}(X))^2}{n}} = \text{RMS of deviations from average}. \nonumber \end{equation}\]
The SD moves avoid the computation problems with the a.a.d.
To calculate the SD, first make this little table, with the list of values, the deviations from the average, and the squares of the deviations.
| \(X\) | \(d\) | \(d^2\) |
|---|---|---|
| 1 | -2 | 4 |
| 2 | -1 | 1 |
| 3 | 0 | 0 |
| 4 | 1 | 1 |
| 5 | 2 | 4 |
Then compute the average of the squares of the deviations, which in this case is 2. Then take the square root of that average, which in this case is about 1.4. Notice that 1.4 is about 1.2 (the a.a.d.). The SD is bounded (weakly) below by the a.a.s., but they’ll usually be close, so we can think of the SD as how far a typical point falls away from the average.
Exercise 5.9 Recall that we define the SD as the RMS of the deviations from average. This definition reveals the steps to compute the SD–just work your way backwards through it. Remember that “R” stands for “root” (or “square root”), “M” stands for “mean” (or “average”), and “S” stands for “square.”
- Find the average.
- Find the deviations from the average.
- Find the S of the deviations from the average.
- Find the MS of the deviations from the average.
- Find the RMS of the deviations from the average. This is the SD.
Solution
- 5,
- 0, -1, 3, -3, 1
- 0, 1, 9, 9, 1
- 4+
- 2
Exercise 5.10 Find the average and SD of the following lists of numbers. Round to two decimal places where necessary.
- 4, 2, 5, 6, 7
- 5, 2, 7, 10, 18
- 87, 123, 47, 62, 153
5.2.4 Other Summaries
The combination of the average and SD is sometimes called a “two-number summary.” This combinations summarizes two important features of a variable using only two numbers. However, there are a few more ways to summarize a variable.
- Minimum or Min: the smallest value of the varaible.
- Maximum or Max: the largest value of the variable.
- 25th Percentile or 1st Quartile: the point which 25% of the data fall below and 75% fall above. To compute, arrange the data from lowest to highest and divide the data into an upper and lower half. If the total number of values is odd, the exclude the median from both halves. The 25th percentile is the median of the lower half.
- 75th Percentile or 3rd Quartile: the point which 75% of the data fall below and 25% fall above. The 75th percentile is the median of the upper half from above.
- Interquartile Range or IQR: the difference between the 75th and 25th percentiles. Intuitively, this is the width of the middle 50% of the data.
We sometimes refer to the median, 25th and 75th percentiles, and minimum and maximum as the five-number summary.
Exercise 5.11 Find the five-number summaries (min, 25th percentile, median, 75th percentile, max) of the following lists of numbers. Also compute the interquartile range. Round to two decimal places where necessary.
- 27, 9, 6, 38, 51, 61, 61, 95, 72
- 33, 6, 0, 28, 21, 5, 86, 67
- 3.57, 1.43, 3.55, 4.82, 2.59, 4.41, 5.14, 2.84, 4.89
- 3.8, 5.9, 6.3, 4.4, 8.3, 9.1, 6.8, 8.1, 5.4, 9.9, 5.8, 7.2
Exercise 5.12 As political scientists, we think a lot about trends in partisanship and voting behavior. Going back to The American Voter, published in 1960, the partisan identity of an individual has seemed to be the dominant factor in determining for whom they cast a vote.
- Before Looking at the figures below, how do you think partisan trends have shifted since the 1950s? Have people become more or less partisan? Why do you think this is? Write down a paragraph or so sketching out your model of partisanship trends.
Exercise 5.13 The figures below are taken from “Partisanship and Voting Behavior, 1952-1996” published in 2000 by Larry Bartels.
Figure 1 shows the proportion of survey respondents who identified as a “Weak Partisan” or a “Strong Partisan”, as well as those who consider themselves “Independent Leaners” or “Pure independents”. An Independent leaner is someone who may not be registered with a party (or does not think of themselves as ‘belonging’ to either party), but nevertheless consistently supports one party over the other.
What proportion of the sample identified as partisan in 1952? Is this more or less than in 1996?
Assume approximately 10,000 people responded to the survey in each year of the survey. Approximately how many more people considered themselves weak identifiers in 1976 than strong identifiers?
Again, assume 10,000 responses in each iteration of the survey: are there more independents in 1996 than partisans? Of the independents, are more “leaners” or “pure”?
Figure 2 shows the proportion of voters and nonvoters who identify themselves with one party or the other.
In general, are those who vote more or less partisan than those who do not?
In terms of partisanship, have voters and nonvoters become more or less similar?
Exercise 5.14 How well do the data above support (or not) your hypothesized model? How would you change your model to accommodate this new information?
- These data end in 1996. If we collected data up to 2020 for figures 1 and 2, do you think the trends you identified would continue? Why?
Some people refer to the “average” as the “mean”. I prefer to avoid this because the “mean” might also refer to the expectation of a random variable. I use “average” and “expected value” to differentiate these two meanings.↩
The median, alternatively, is not easy to compute and quite difficult to work with mathematically.↩


Power Usage Effectiveness
Power usage effectiveness (PUE) is a metric developed by The Green Grid that measures data center efficiency. It’s a ratio of the total power consumed by a data center to the energy delivered to the servers. The ideal PUE is 1.0, meaning 100% of the data center energy powers the IT equipment. Best practice PUE for the IT sector is 1.5.