search this blog

Thursday, June 9, 2016

The discrepancy


I posted this in the comments in the previous blog entry [here], but it hasn't received the attention that it deserves, so it's now getting an entry all of its own. I'm also posting Matt's reply, since he motivated me to look at the formal stats from Hofmanova et al. in more detail. I'd like to get to the bottom of this. Any ideas?

Davidski said...

Matt,

How would you interpret these sets of f4 and D statistics?

The f4 stats are from Hofmanova et al., while the D stats were run by me. The first set of D stats uses the highest quality Anatolia Neolithic sample from Barcin from Mathieson et al. and CHG genotypes from Fu Q et al., and the second uses the same Barcin sample plus CHG genotypes from Jones et al.

Also, keep in mind that, as far as I can tell, the Barcin genomes from Hofmanova et al. and Mathieson et al. date to the same period.

f4 Corded_Ware_LN Bar8 Satsurblia Khomani -0.0367 -8.145
f4 Corded_Ware_LN Bar8 Kotias Khomani -0.0193 -3.437

f4 Spain_MN Bar8 Satsurblia Khomani -0.0327 -5.385
f4 Spain_MN Bar8 Kotias Khomani -0.0182 -3.136

versus...

D Corded_Ware_Germany BAR20_I0709 Satsurblia Khomani 0.0215 4.299
D Corded_Ware_Germany BAR20_I0709 Kotias Khomani 0.0205 4.408

D Iberia_MN BAR20_I0709 Satsurblia Khomani -0.0003 -0.06
D Iberia_MN BAR20_I0709 Kotias Khomani -0.0017 -0.339

...

D Corded_Ware_Germany BAR20_I0709 Satsurblia2 Khomani 0.0224 4.38
D Corded_Ware_Germany BAR20_I0709 Kotias2 Khomani 0.0226 5.168

D Iberia_MN BAR20_I0709 Satsurblia2 Khomani 0.0068 1.222
D Iberia_MN BAR20_I0709 Kotias2 Khomani 0 0.004

Clearly, something's horribly wrong. If I made a mistake, my apologies. But I'm pretty sure I didn't make any mistakes. I checked the datasets that I'm using for consistency with the f4 and D stats published in Mathieson et al. and Fu et al., so I can say with confidence that my D stats should not be much different from correctly run f4 stats using the same ancient samples.

June 8, 2016 at 7:06 PM

Matt said...

@ Davidski, yeah I see what's going on there with the D stats giving a result we would expect from previous work - Anatolia Neolithic and Iberia_MN equally related to CHG, Corded Ware more to CHG - with the stats from this paper being different - Bar8 being more related to CHG than Iberia_MN, and Bar8 even more strongly related to CHG relative to Corded_Ware, also implying Iberia_MN more related to CHG than Corded Ware is. I don't know that there's anything about f4 vs D stats themselves that would explain that difference, and as you say, yours are consistent with the previously published.

This is really stuff that should have been in and the resolved in the early print. That's the whole point of the process!

June 9, 2016 at 12:13 AM

Update 12/06/2016: Here are f4 stats using the same data as for the D stats above. They look basically the same as the D stats.

Corded_Ware_Germany BAR20_I0709 Satsurblia Khomani 0.002073 4.294
Corded_Ware_Germany BAR20_I0709 Kotias Khomani 0.001997 4.402

Iberia_MN BAR20_I0709 Satsurblia Khomani -0.00003 -0.06
Iberia_MN BAR20_I0709 Kotias Khomani -0.000157 -0.34

Corded_Ware_Germany BAR20_I0709 Satsurblia2 Khomani 0.0021 4.388
Corded_Ware_Germany BAR20_I0709 Kotias2 Khomani 0.002031 4.971

Iberia_MN BAR20_I0709 Satsurblia2 Khomani 0.000613 1.221
Iberia_MN BAR20_I0709 Kotias2 Khomani 0.000002 0.004

20 comments:

Samuel Andrews said...

Have you tested the accuracy of any of their other D-stats. It'd be helpful to either do that or wait for the genomes from this paper to be available and then do analysis of them. Without that we can't say anything about these new genomes besides that they were EEF. I've heard Kum6 from circa 3000 BC had lots of CHG affinity, is this true?

Davidski said...

Have you tested the accuracy of any of their other D-stats.

They didn't run any D-stats. They ran f4-stats using their own program.

But yeah, I did test a whole bunch of their f4-stats that didn't involve the CHG genomes, and these matched my D-stats very closely.

Alberto said...

Yes, there is something strange in the stats. But I see that in the Admixture runs, both Bar8 and Bar31 (as well as the samples from Greece) look quite distinct from EEFs and quite CHG-like. So it doesn't seem to be a problem specific to F4 stats.

I hope the samples are available soon. It's going to be interesting to see how they behave in more tests and how they can expand our understanding of Neolithic Europe.

Karl_K said...

It is published, so you should be able to ask for the data. A little analysis will sort this out in no time. Sample names get switched all the time.

Roy King said...

@Davidski,
I know you tend to anathemize semantics and philosophical analysis, but I just wanted to add my two cents about scientific methodology/culture based on years of reviewing papers, grants and investigating scientific misconduct as a dean.
In order to refute the f4 results of the paper, you should first use their data directly and replicate their f4 analyses, including the two Kumtepe samples and the CHG samples (or fail to replicate). If you can't replicate, you should write a letter to PNAS asking about either their retracting their paper or having PNAS publish your non-confirming analyses. Difference in semantic interpretation of particularly minor theses from the paper, like the predicting of modern DNA from ancient samples, although revealing, do not justify retraction, but necessitate explication. By the way, we all get PHDs, Doctors of Philosophy, as a homage to the philosophical roots of our scholarship.

Nirjhar007 said...

I think that David was too emotional .Something made him very upset.

Davidski said...

Roy,

The results and samples are being checked by some very capable people. Not sure how long it'll take, but probably not too long.

Nirjhar,

I get very angry and aggressive when I see major journals publishing pop gen stuff that should never be published.

Davidski said...

Here are some f4 stats:

f4 Corded_Ware_Germany Anatolia_Neolithic CHG Chimp 0.002396 9.226 574503
f4 Corded_Ware_Germany Anatolia_Neolithic CHG Mbuti 0.002084 9.113 587219
f4 Yamnaya_Samara Anatolia_Neolithic CHG Chimp 0.002836 10.59 574578
f4 Yamnaya_Samara Anatolia_Neolithic CHG Mbuti 0.002661 10.793 587274
f4 Esperstedt_MN Anatolia_Neolithic CHG Chimp 0.00044 0.98 494134
f4 Esperstedt_MN Anatolia_Neolithic CHG Mbuti -0.000031 -0.083 504920
f4 Baalberge_MN Anatolia_Neolithic CHG Chimp 0.001039 2.326 286950
f4 Baalberge_MN Anatolia_Neolithic CHG Mbuti 0.0003 0.856 293152

bau said...

Davidski,
where did you read Hofmanova et al used their own program for f4? They claim that used ADMIXTOOLS, or did I miss something?

"f3 and f4 statistics and the associated Z-scores (via block jackknife with default options) were calculated using the ADMIXTOOLS package"

Davidski said...

Someone told me, but never mind, if it was ADMIXTOOLS then OK.

bau said...

But, sorry if I ask, Is it possible to estimate f4 with admixtools? I didn't fine the software in the readme provided with the suite, is it something enclosed in qpf4ratio?

Davidski said...

Not in the current version it's not. Next version maybe, so maybe they had the new version, but I doubt it.

Nick Patterson (Broad) said...

@bau

Sure you can get f4-stats; use qpDstat and code
f4mode: YES ; this has been a feature for a long time...

By the way, expect a new release of ADMIXTOOLS very soon

(probably this coming week) which will have a public domain version of
qpGraph. Check the Reich Lab software page.

bau said...

@Nick
Thank you very much for this! I didn't find it before.

bau said...

@nick
Actually i couldn't find the f4mode oprtion in the readme so it's not clear if I should edit the par file or the source code and recompile. Could you please help me?

Davidski said...

Put this in the parfile...

f4mode: YES

Davidski said...

Used to do f4 with fourpop. This makes it easier.

f4 stats for the same data as for the D stats above...

Corded_Ware_Germany BAR20_I0709 Satsurblia Khomani 0.002073 4.294
Corded_Ware_Germany BAR20_I0709 Kotias Khomani 0.001997 4.402

Iberia_MN BAR20_I0709 Satsurblia Khomani -0.00003 -0.06
Iberia_MN BAR20_I0709 Kotias Khomani -0.000157 -0.34

Corded_Ware_Germany BAR20_I0709 Satsurblia2 Khomani 0.0021 4.388
Corded_Ware_Germany BAR20_I0709 Kotias2 Khomani 0.002031 4.971

Iberia_MN BAR20_I0709 Satsurblia2 Khomani 0.000613 1.221
Iberia_MN BAR20_I0709 Kotias2 Khomani 0.000002 0.004

Tobus said...

What's the difference between f4 and D-stat? I thought they were essentially the same thing...?

Davidski said...

They're similar, but not the same. This is from Nick's paper from a few years ago...

As mentioned earlier, D-statistics are very similar to the 4-population test statistics introduced in REICH et al. (2009). The primary difference is in the computation of the denominator of D. For statistical estimation, and testing for ‘treeness’, the D-statistics are preferable, as the denominator of D, the total number of ‘ABBA’ and ‘BABA’ events, is uninformative for whether a tree phylogeny is supported by the data, while D has a natural interpretation: the extent of the deviation on a normalized scale from -1 to 1.

http://www.genetics.org/content/early/2012/09/06/genetics.112.145037

I always thought D stats were better. Also, running f4 with fourpop was a bit of a pain in the butt.

Davidski said...

An updated version of ADMIXTOOLS, version 4.1, has just appeared at GitHub.

https://github.com/DReichLab/AdmixTools