This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Calculate pairwise Fst from Plink (.bim .bed .fam) files in R?

Hi there. I'm new to bioinformatics analysis of SNP data. I downloaded a SNP dataset generated in Plink and I have successfully loaded it into R as a SnpMatrix so that I can calculate pairwise Fst using the snpStats package:

geno<-read.plink("dataset1b.bed", "dataset1b.bim", "dataset1b.fam", na.strings = c("0", "-9"), sep = "." , select.subjects = NULL, select.snps = NULL)
head(geno)

$genotypes
A SnpMatrix with  67 rows and  96899 columns
Row names:  1 ... 67 
Col names:  1_4683 ... 31_24666833 

$fam
   pedigree member father mother sex affected
1         1      1     NA     NA  NA       NA
2         2      1     NA     NA  NA       NA
3         3      1     NA     NA  NA       NA
4         4      1     NA     NA  NA       NA
5         5      1     NA     NA  NA       NA
6         6      1     NA     NA  NA       NA
7         7      1     NA     NA  NA       NA .... 

$map
          chromosome  snp.name cM position allele.1 allele.2
1_4683             1    1_4683 NA     4683        C        T
1_4739             1    1_4739 NA     4739        A        G
1_4793             1    1_4793 NA     4793        C        T
1_4805             1    1_4805 NA     4805        T        C
1_160088           1  1_160088 NA   160088        A        G
1_160094           1  1_160094 NA   160094        A        G
1_160102           1  1_160102 NA   160102        A        G
1_160106           1  1_160106 NA   160106        T        C
1_168229           1  1_168229 NA   168229        T        C
1_168262           1  1_168262 NA   168262        T        C

So now I need to use Fst(snps, group, pairwise=TRUE) Arguments: snps: an object of class SnpMatrix or XSnpMatrix containing the SNP data group: a factor (or object than can be coerced into a factor), of length equal to the number of rows of snps, giving the grouping or rows for which the Fst is to be calculated pairwise: if TRUE, the within-group variances are weighted according to the number of possible within-group pairwise comparisons of chromosomes. If FALSE, the default value, weights are simply the number of chromosomes in each group.

However, I am not sure what to use for "group". The fam file contains info on the samples, but when I do

f1<-Fst(geno, geno$fam, pairwise=T)

I get

Error in Fst(geno, geno$fam, pairwise = T) : 
  STRING_ELT() can only be applied to a 'character vector', not a 'NULL'

Any thoughts? I see that "group" says that the length must be equal to the number of rows of SNPs, but the Plink output is all I have to work with.

Thanks in advance!

plink r snp snpstats fst

0 answers

No answers yet.

Log in to answer this question.