PCAngsdTutorial: Difference between revisions
Line 10: | Line 10: | ||
<code> | <code> | ||
PCAngsd=~/Software/ | PCAngsd=~/Software/pcangsd | ||
</code> | </code> | ||
Revision as of 23:07, 29 July 2019
We will go through a simple and more complex example on how to use PCAngsd with visualization of the data.
PCA with admixture aware priors
This example will perform a PCA analysis on 1000 genotype likelihoods.
Every time you open a new terminal window, set the paths to the program and the input file.
Set the path to PCAngsd
PCAngsd=~/Software/pcangsd
Test the link
ls $PCAngsd
Create directories
Create the directories that will be used for working:
mkdir Demo
cd Demo
mkdir Data
mkdir Results
Set the paths to your local directories
IN_DIR=Demo/Data
OUT_DIR=Demo/Results
Test the links
ls $IN_DIR
ls $OUT_DIR
Download the beagle genotype likelihood input file
PCAngsd uses Genotype Likelihoods (GLs) in .beagle format as input. The input file has been created for you.
Download the files and move them to your input folder (for example, $IN_DIR):
*ANDERS*
wget popgen.dk/software/download/NGSadmix/data/Demo1input.gz
wget popgen.dk/software/download/NGSadmix/data/Demo1pop.info
mv Demo1input.gz $IN_DIR
mv Demo1pop.info $IN_DIR
*ANDERS*
View the genotype likelihood beagle 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/Demo1input.gz | head -n 10 | cut -f 1-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/Demo1input.gz | wc -l
View population information file
To view a summary of the population information file, cut the first column, sort and count:
cut -f 1 -d " " $IN_DIR/Demo1pop.info | sort | uniq -c
Lets make a population label file and place it in the output directory
cut -f1 -d" " $IN_DIR/Demo1pop.info > $OUT_DIR/poplabel
Load the Python module
PCANGSD='python ~/software/pcangsd/pcangsd.py'