NgsAdmixTutorial: Difference between revisions
No edit summary |
|||
Line 3: | Line 3: | ||
==Example of NGSadmix - Simple== | ==Example of NGSadmix - Simple== | ||
In our first example, we will | In our first example, we will infer admixture proportions for low depth NGS data using a small dataset from 30 human samples. | ||
Set directories to all required programs and the data you will use depending on where you installed them every time you open a new terminal window | # Set directories to all required programs and the data you will use depending on where you installed them every time you open a new terminal window | ||
# Set path to NGSadmix | # Set path to NGSadmix | ||
Line 14: | Line 14: | ||
# Set path to the directory with the input files | # Set path to the directory with the input files | ||
AA=/../home/albrecht/PhDCourse | AA=/../home/albrecht/PhDCourse | ||
#test the link | #test the link | ||
ls $AA | ls $AA | ||
Line 19: | Line 20: | ||
https://github.com/mfumagalli/ngsTools/blob/master/Scripts/data.sh | https://github.com/mfumagalli/ngsTools/blob/master/Scripts/data.sh | ||
Create the directories that will be used for working: | #Create the directories that will be used for working: | ||
mkdir Demo | mkdir Demo | ||
cd Demo | cd Demo | ||
Line 33: | Line 34: | ||
ls $OUT_DIR | ls $OUT_DIR | ||
NGSadmix uses Genotype Likelihoods (GLs) in .beagle format as input. The input file has been created for you. | #Input file | ||
NGSadmix uses Genotype Likelihoods (GLs) in .beagle format as input. The input file has been created for you. | |||
We will use a very reduced data set: | |||
-10 individuals from each population: 10 from Nigeria (YRI), 10 from Japan (JPT) and 10 with European Ancestry (CEU). | |||
-a very reduced genome 30 x 100k random regions across the autosomes | |||
-each individual is sequenced at 2-6X | |||
<table class="muse-table" border="2" cellpadding="5"> | <table class="muse-table" border="2" cellpadding="5"> | ||
Line 52: | Line 56: | ||
</table> | </table> | ||
#Population information file. | |||
A file with labels that indicate which population the individuals are sampled from is also provided as the population information. | |||
#Copy the files to your input folder ($IN_DIR): | |||
cp popgen.dk/software/download/NGSadmix/data/Demo1GLs $IN_DIR | |||
cp popgen.dk/software/download/NGSadmix/data/pop.info $IN_DIR | |||
We need to prepare the population information file by cutting the first column, sorting and counting: | #We need to prepare the population information file by cutting the first column, sorting and counting: | ||
cut -f 1 -d " " $IN_DIR/pop.info | sort | uniq -c | cut -f 1 -d " " $IN_DIR/pop.info | sort | uniq -c | ||
Revision as of 22:09, 11 July 2019
We will go through a simple and more complex example on how to use NGSadmix with visualization of the data.
Example of NGSadmix - Simple
In our first example, we will infer admixture proportions for low depth NGS data using a small dataset from 30 human samples.
- Set directories to all required programs and the data you will use depending on where you installed them every time you open a new terminal window
- Set path to NGSadmix
NGSADMIX=~/Software/NGSadmix
- test the link
ls $NGSADMIX
- Set path to the directory with the input files
AA=/../home/albrecht/PhDCourse
- test the link
ls $AA To create the data files, please go the the following link: https://github.com/mfumagalli/ngsTools/blob/master/Scripts/data.sh
- Create the directories that will be used for working:
mkdir Demo cd Demo mkdir Data mkdir Results
- Set paths to local directories
IN_DIR=~/Demo/Data OUT_DIR=~/Demo/Results
#Test the links ls $IN_DIR ls $OUT_DIR
- Input file
NGSadmix uses Genotype Likelihoods (GLs) in .beagle format as input. The input file has been created for you. We will use a very reduced data set: -10 individuals from each population: 10 from Nigeria (YRI), 10 from Japan (JPT) and 10 with European Ancestry (CEU). -a very reduced genome 30 x 100k random regions across the autosomes -each individual is sequenced at 2-6X
CEU | Europeans (mostly of British ancestry) |
JPT | East Asian - Japanese individuals |
YRI | West African - Nigerian Yoruba individuals |
- Population information file.
A file with labels that indicate which population the individuals are sampled from is also provided as the population information.
- Copy the files to your input folder ($IN_DIR):
cp popgen.dk/software/download/NGSadmix/data/Demo1GLs $IN_DIR cp popgen.dk/software/download/NGSadmix/data/pop.info $IN_DIR
- We need to prepare the population information file by cutting the first column, sorting and counting:
cut -f 1 -d " " $IN_DIR/pop.info | sort | uniq -c
Now let’s analyze the input file:
In general, the first three columns of a beagle file contain marker name and the two alleles, allele1 and allele2, present in the locus (in beagle A=0, C=1, G=2, T=3). All following columns contain genotype likelihoods (three columns for each individual: first GL for homozygote for allele1, then GL for heterozygote and then GL for homozygote for allele2). Note that the GL values sum to one per site for each individual. This is just a normalization of the genotype likelihoods in order to avoid underflow problems in the beagle software, but it does not mean that they are genotype probabilities.
In order to see the first 10 columns and 10 lines of the input file, type: gunzip -c $IN_DIR/Demo1GLs | head -n 10 | cut -fl-10 | column -t
Use this command to count the number of lines of the input file. The number of lines, indicates the number of loci for which there are GLs plus one (as the command includes the count of the header line): gunzip -c $IN_DIR/Demo1GLs | wc -l
To run an analysis of the GLs with NGSadmix, assuming the number of ancestral populations is K=3, type the following command:
$NGSadmix -likes $IN_DIR/Demo1GLs.beagle.gz -K 3 -minMaf 0.05 -seed 1 -o $OUT_DIR/Demo1NGSadmix
The output:
The analysis performed by NGSadmix produces 4 files:
Log likelihood of the estimates A .log file that summarizes the run. Let’s take a look at the log file to determine the log likelihood of the estimates achieved by NGSadmix which is called the “best like” in the file: cat Demo1NGSadmix.log
Estimated allele frequency A zipped .fopt file, that contains an estimate of the allele frequency in each of the 3 assumed ancestral populations (there is a line for each locus). We can use this file to obtain the estimated allele frequency of the first 5 SNPs (one per line) of the three assumed ancestral populations, by typing the following command: zcat Demo1NGSadmix.fopt.gz | head -n 5
Estimated admixture proportions A .qopt file, that contains an estimate of the individual's ancestry proportion from each of the three assumed ancestral populations (there is a line for each individual). To obtain the estimated admixture proportions for the first 5 individuals, type the following command: head -n 5 Demo1NGSadmix.qopt
DemoNGSadmix.filter: if the filter was used, it will show the sites that were left out. If the filter is used, it will show the sites that were left out or in… CHECK
Follow these instructions to make a simple plot of the estimated admixture proportions for all individuals in R:
Type “R” in the terminal and press enter. Paste the following code into R:
- set up the working directories
- Fill up a table with the IDs of the population information for each individual
pop<-read.table("pop.info", as.is=T)
- Read inferred admixture proportions
q<-read.table("Demo1NGSadmix.qopt")
- Plot them (ordered by population)
ord = order(pop) par(mar=c(7,4,1,1)) barplot(t(q)[,ord],col=c(2,1,3),names=pop[ord],las=2,ylab="Admixture proportions",cex.names=0.75)
The y-axis of the plot show the admixture proportions and the individuals in the samples are plotted in the x-axis. Each color represents a different ancestral population. The proportion of each color shows the different admixture of the individuals for each ancestral population. The plot is sorted by the population of origin of each individual in the sample, and therefore, it shows blocks with prevalence of a certain color, which represents the population to which each individual belongs.
NB As you could tell from the number of loci included in the analysis, the above analysis is based on data from very few loci (actually we on purpose only analyzed data from a small part of the genome to make sure the analysis ran fast). Results from an analysis of data from the entire genome can be seen here.