Plink --score issues

[Bs He]

The objective is to calculate the genetic risk score, given the genotyping data, effective allele and effective size. Totally 1.5 million SNPs are included in the bfile.bim, and they all have variant ID with form Chr:BP; e.g., 1:12345. However, after I submit the code

/projects/bsi/gentools/bin/plink2 --bfile GenotypingData --score ScoreFile header sum --threads 6

Then I have the result:

     FID        IID  PHENO    CNT      CNT2    SCORESUM
       01        01     -9    3067556 2438692    9.411
       02        02     -9    3067556 2440321  9.16466
       03        03     -9    3067556 2440784  9.50342
       04        04     -9    3067556 2443276   10.615

About 1.5 million snps are concerned in the score. So I can understand thant CNT is approximately 2*1.5 million (for diploid genomes). However, how to understand CNT2? What does "Sum of named allele counts" mean as in plink.profile description.

Second question: in pline allelic scoring function, it says "Also, note that scores are multiplied by 0..1 dosages, not 0..2 diploid allele counts, unless the 'double-dosage' modifier is present". What does this mean? What is the difference between 0..1 dosages and 0..2 diploid allele counts?

Third question: To use --score ScoreFile function, in ScoreFile, we need SNP ID, effect allele and effect size. Say, if the effect allele is the minor allele, then the score contributed by the SNP should be effective size * dosage (effect allele); however, if the effect allele is not the minor allele, then how to calculate the contribution of the SNP. It should be effect size * (2 - dosage (effect allele)), or effect size * - dosage(effect allele)?

I think the third question is also related with the second one. Hope somebody can help.

[Christopher Chang]

1. Each ScoreFile line names a variant, and only one of its two alleles.  CNT2 is the sum of the counts of only the named alleles.

2. This is irrelevant if you aren't using --score with --dosage.  (There is no longer any reason to do so; use plink 2.0 --import-dosage instead for that use case.)

3. The same formula is used regardless of whether the effect allele is major or minor.  Non-effect alleles have no impact on score unless you include the 'center' modifier (or plink 2.0 --score's 'variance-standardize' modifier).

[Bs He]

Thanks for such reply. However, I am not sure I understand what you said.

1. Each ScoreFile line names a variant, and only one of its two alleles.  CNT2 is the sum of the counts of only the named alleles. 

Just a quick question, how does Plink get CNT2=2438692 this result for the 1st subject. I mentioned "Totally 1.5 million SNPs are included in the bfile.bim, and they all have variant ID with form Chr:BP; e.g., 1:12345. ", which means the bfiles GenotypingData.* were generated by using --extract 1.5millionSNPs_list. Therefore, the bfile contains and just contains the exactly the same SNPs used in the score.

I checked the GenotypingData.raw file from using --recode the bfiles. How should I get 2438692 by operating on the rows?  

3. The same formula is used regardless of whether the effect allele is major or minor.  Non-effect alleles have no impact on score unless you include the 'center' modifier (or plink 2.0 --score's 'variance-standardize' modifier).

Do you mean we alway use effective size * dosage (effect allele)? Say, if the effect allele is A, but the minor allele is G. Since we only have the dosage of G, then I think the dosage of A should be (2-dosage(G)), right? Then the score should be beta*(2-dosage(G))?

Sorry to disturb you. I am a complete novice to genetics, and I really appreciate your help.

[Bs He]

For the first question, what makes the inconsistency between #(SNPs in score) and CNT2? Or 2*#(SNPs in score) and CNT2?

Thanks.

On Tuesday, August 28, 2018 at 3:08:16 PM UTC-5, Christopher Chang wrote:

[Christopher Chang]

Please reread my previous response.  The SCORE FILE only names ONE of the two alleles for each variant.

[Bs He]

So I think you actually mean is, say for one variant from the bim file

1 1:752721 0 752721 A G

So if the effect allele is the ref allele (i.e., G here), then the corresponding count in CNT2 is dosage(G); otherwise, if the effect allele is the alternative allele (i.e., A here), then the corresponding count in CNT2 is dosage(A), right?

Why I made such guess? Because I calculated the count of ref allele for each subject, and I found such count is very close to CNT2 (slight less than CNT2), and meanwhile most effect alleles are the ref alleles. Am I correct? 

[Christopher Chang]

You don't seem to understand the words "SCORE FILE", so the communication barrier is too high for me to help you in any reasonable amount of time.  Sorry.

[Bs He]

Well... actually my guess in the last post is correct, and I just verified it form my code. For the previous miscommunication, I should have clearly in each post indicated that THE ALLELE NAME IN MY BIFILE IS EXACTLY THAT IN THE SCOREFILE. 

For people who may have the same silly problem as me, I straightforwardly give the answer here.

Suppose your ScoreFile has 1 line, 

rsID               effect_allele          beta

1:2245570        G            -0.0276009 

your GenotypingData.bim has a corresponding line ( ignore the ambiguous strands problem here)

1   1:2245570   2245570   G   C 

Therefore the ref allele (major allele) is the effect allele.

Suppose you have a patient with FID 16214852, in the result PLINK.PROFILE file, you find the result

   16214852   16214852     -9      2      0        0

CNT2 = 0 AND the score = 0. 

Now we try to figure out how Plink get the result. Let's recode the BFILES to .raw file, and we extract the dosage of patient 16214852 and the corresponding dosage of 1:2245570. I use awk to do that, and you can use R or cpp. The result is 

1:2245570_C

2

Therefore, the dosage of the effect allele G = 0, that is why CNT2 = 0. So in human easy-understanding language, CNT2 = SUM OF DOSAGE(EFFECT ALLELE).

To be honest, I would like to know why in Plink help of Allelic scoring, the term of "effect allele" is not explicitly mentioned. Hope can be more friendly to very beginning users.

[Christopher Chang]

--score has other uses, such as PCA projection.  "Effect allele" is inaccurate there.  "Allele named in the --score file" is always accurate.

[Bs He]

Nice, thanks.