PDA

View Full Version : Combining different Reich dataset samples in a D-stat run



Bas
10-13-2017, 09:44 PM
I've been able to fun a few D-stats with the Reich datasets but would like to include one or two samples from different Reich sets to run. The Reich datasets are of course all .geno/.ind/.snp. Is this possible?

anglesqueville
10-14-2017, 07:39 AM
I've been able to fun a few D-stats with the Reich datasets but would like to include one or two samples from different Reich sets to run. The Reich datasets are of course all .geno/.ind/.snp. Is this possible?

Haha! That is the fatal slope :). You have now to merge datas. Inside Admixtools there is a program "mergeit" that makes that, but it is imo very unreliable. You have to install Plink, and work with bed/bim/fam files. For Plink learning there are many internet sites, but this one is excellent: https://www.cog-genomics.org/plink2/data. To convert plink format to eigenstrat you'll use "convertf" from Admixtools. Try all this, and ask if you experience problems. I have worked with many files these last months, and I'm now a little aware of those problems. But of course I'm still a modest amateur, and there are actual pros here, who will be able to help you far better than I ( Kurd, Generallissimo, etc).

lukaszM
10-14-2017, 05:51 PM
BTW is it possible to convert geno/.ind/.snp. to plink format (bim,bed,fam)?

Chad Rohlfsen
10-14-2017, 11:05 PM
You use convertf from AdmixTools.

Bas
10-18-2017, 05:10 PM
Haha! That is the fatal slope :). You have now to merge datas. Inside Admixtools there is a program "mergeit" that makes that, but it is imo very unreliable. You have to install Plink, and work with bed/bim/fam files. For Plink learning there are many internet sites, but this one is excellent: https://www.cog-genomics.org/plink2/data. To convert plink format to eigenstrat you'll use "convertf" from Admixtools. Try all this, and ask if you experience problems. I have worked with many files these last months, and I'm now a little aware of those problems. But of course I'm still a modest amateur, and there are actual pros here, who will be able to help you far better than I ( Kurd, Generallissimo, etc).

Yes,haha! I'm having a few problems already. I tried using mergeit and I don't know if it is the folder I am putting it in (downloaded EIG-7.2.1 and copied the files I want to merge inside and then that whole folder inside AdmixTools-master-my run path is something like /home/liveuser/AdmixTools-master/EIG-7.2.1/bin/mergeit -p parfile...) but I keep getting the following error message:


error while loading shared libraries: libgfortran.so.3: cannot open shared object file: No such file or directory

Weird because I thought that the program would have everything necessary and a shame because I already have all the necessary files in .geno/.ind/.snp format >:(

anglesqueville
10-18-2017, 08:00 PM
Likely you don't have fortran3 installed. Try to install libgfortran3:i386 ( yum install libgfortran3:i386 ). That said I'm not a fan of mergeit. Install Plink and use --bmerge ( doc on cog-genomics). You'll likely experience problems ( very likely snps with 3+ alleles for example). If you are lost, call for help.

Bas
10-22-2017, 08:40 PM
Likely you don't have fortran3 installed. Try to install libgfortran3:i386 ( yum install libgfortran3:i386 ). That said I'm not a fan of mergeit. Install Plink and use --bmerge ( doc on cog-genomics). You'll likely experience problems ( very likely snps with 3+ alleles for example). If you are lost, call for help.

Thanks, finally got time and managed to convert two datasets to ped/map and then in PLINK to bed/bim/fam, then merge them both like this. No problems with that or converting the merged file back to ped/map but when converting back to EIGENSTRAT in convertf I got the message


fatalx:
duplicate sample: I0011
Aborted (core dumped)

What should I do here?

anglesqueville
10-22-2017, 08:43 PM
Thanks, finally got time and managed to convert two datasets to ped/map and then in PLINK to bed/bim/fam, then merge them both like this. No problems with that or converting the merged file back to ped/map but when converting back to EIGENSTRAT in convertf I got the message



What should I do here?

Duplicate sample? Never experienced this problem. If I'm well understanding, you have 2 I0011 in your .fam?

Bas
10-23-2017, 10:05 PM
Duplicate sample? Never experienced this problem. If I'm well understanding, you have 2 I0011 in your .fam?

Yes, that was the problem...along with 72 other duplicate entries,haha. That got sorted but then I ran into another problem-the resulting EIG files weren't working.

Normally the .GENO files from Reich cannot be opened in my Fedora. No problems with using them however but when I right click for the properties on that, the file shows as 'program (application/octet-stream)' type. But when I changed my merged PED/MAP files back to EIGENSTRAT in convertf, the resulting .GENO file was just a plaintext file, which I could open and the whole text file full of 2099229292909 etc,across the whole page . When I went to run dstats, the program hung for about 20 minutes before giving me a 'killed' message.

I get the feeling due to the difference between my new .GENO file and the original Reich .GENO that the problem is something to do with this.

Between converting original EIG datasets in convertf and then PLINK, merging and then converting back to EIG I didn't do anything extra to the files (like marker prune). I don't know if I should be doing something before I merge or after I merge or before I convert back to EIG? Or if this is even the cause of the problem?

anglesqueville
10-24-2017, 06:15 AM
It's difficult to remotely solve such problems. Perhaps if you tell me which datas you wanted to merge I'll can try to do it myself and see what happens.

Bas
10-24-2017, 02:25 PM
It's difficult to remotely solve such problems. Perhaps if you tell me which datas you wanted to merge I'll can try to do it myself and see what happens.

That would be extremely kind of you, it was the Minoan/Mycenaean dataset from Lazaridis 2017 with the Haak 2015 fully public dataset.

anglesqueville
10-24-2017, 03:42 PM
Fortunately I have those two datas, converted to EIG, on my disk. Many samples are in both. Well... here is what I would do, but there is likely other solutions ( and perhaps better ones). Assuming that Haak is full.bed/bim/fam, and the other MinMyc.bed/bim/fam:
1) With the text editor, extract from full.fam the samples which are in MinMyc.fam. Call the file "common_samples"
2) In Plink: ./plink --bfile full --remove common_samples --make-bed --allow-no-sex --out full_1
3) Now you can merge full_1 and MinMyc. Personally ( again I'm not claiming it's the best way) I use to:
4) If the fams have samples numbers in the first column, I modify the numerotation of the smallest fam in order to get a logical global numerotation: 1....x, x+1....y.
5) Using the flag "--extract" I work in order to get two plink files with the same SNPs pannel. For example, assuming that there are more SNPs in MinMyc ( I believe it's the case), I make:
./plink --bfile MinMyc --write-snplist ----------------> you get a file "plink.snplist"
./plink --bfile full_1 --extract plink.snplist --make-bed --allow-no-sex --out full_2
6) You can now merge full_2 and MinMyc:
./plink --bfile minMyc --bmerge full_2.bed full_2.bim full_2.fam --make-bed --allow-no-sex --out new

If you are lucky, you'll get new.bed/bim/fam. If you are not ( the most probable, according to my personnal experience), you will experience multiallelic snps ( and get a file called something like "new-merge.missnp"). In that case, you can try to flip the bad snps, or roughly exclude them:
7)./plink --bfile full_2 --exclude new-merge.missnp --make-bed --allow-no-sex --out full_3
8) ./plink --bfile minMyc --bmerge full_3.bed full_3.bim full_3.fam --make-bed --allow-no-sex --out new

Now you have likely your new.bed/bim/fam, and you can convert it in EIG

Bas
10-25-2017, 03:15 PM
Fortunately I have those two datas, converted to EIG, on my disk. Many samples are in both. Well... here is what I would do, but there is likely other solutions ( and perhaps better ones). Assuming that Haak is full.bed/bim/fam, and the other MinMyc.bed/bim/fam:
1) With the text editor, extract from full.fam the samples which are in MinMyc.fam. Call the file "common_samples"
2) In Plink: ./plink --bfile full --remove common_samples --make-bed --allow-no-sex --out full_1
3) Now you can merge full_1 and MinMyc. Personally ( again I'm not claiming it's the best way) I use to:
4) If the fams have samples numbers in the first column, I modify the numerotation of the smallest fam in order to get a logical global numerotation: 1....x, x+1....y.
5) Using the flag "--extract" I work in order to get two plink files with the same SNPs pannel. For example, assuming that there are more SNPs in MinMyc ( I believe it's the case), I make:
./plink --bfile MinMyc --write-snplist ----------------> you get a file "plink.snplist"
./plink --bfile full_1 --extract plink.snplist --make-bed --allow-no-sex --out full_2
6) You can now merge full_2 and MinMyc:
./plink --bfile minMyc --bmerge full_2.bed full_2.bim full_2.fam --make-bed --allow-no-sex --out new

If you are lucky, you'll get new.bed/bim/fam. If you are not ( the most probable, according to my personnal experience), you will experience multiallelic snps ( and get a file called something like "new-merge.missnp"). In that case, you can try to flip the bad snps, or roughly exclude them:
7)./plink --bfile full_2 --exclude new-merge.missnp --make-bed --allow-no-sex --out full_3
8) ./plink --bfile minMyc --bmerge full_3.bed full_3.bim full_3.fam --make-bed --allow-no-sex --out new

Now you have likely your new.bed/bim/fam, and you can convert it in EIG

Cheers! Going to try this. Meanwhile, the problem I mentioned seems to have resolved itself. I think it was just a memory issue on my computer because of the large dataset. But that one I merged and the new Tollense data are running without problems so far. But I'm going to try to get cleaner runs with the way you mentioned...and then eventually learn other stuff like qpAdm B)

anglesqueville
10-25-2017, 03:35 PM
I think it was just a memory issue on my computer because of the large dataset
Possible, in particular if you have installed Linux in a virtual computer ( it's my case, I use Oracle VirtualBox). You'll have likely to allocate a rather large amount of RAM to Linux. I experienced myself many memory issues, in particular working with vcf.

But I'm going to try to get cleaner runs with the way you mentioned
Yes, I'm a maniac of clean files, in particular when I merge. But I'm perhaps wrong, and it's perfectly possible that some steps of the procedures I use are needless. After all, I'm myself an amateur, and an autodidact. Anyway, I hope that I've been useful for you, because it's always an excellent thing when someone tries to learn new tools.:thumb:

Bas
10-25-2017, 03:57 PM
Yes, I've actually been running from USB, it does some strange things like files suddenly turning read-only if you manipulate the files from the home directory but once I avoided that by running straight from USB instead of copying files, it's fine. I think the files from the European Nucleotide are VCF aren't they.

You've helped me so much and it's great to have people like you willing to take the time to show others the ropes. :beerchug: