Discrepancy between negative and positive LMM effect conversions

[kheilbron]

Hi Luke, thanks for an excellent paper and repo! I wanted to flag some odd behaviour I've been seeing with the function for converting LMM effects to odds ratios. Flipping the sign of the LMM effect should flip the sign of the resulting logOR, but the absolute value should be unchanged. However, I've found that positive and negative LMM effects get converted differently. Here is a minimal example for you:

df <- data.frame( BETA = c( -0.0075, 0.0075 ), 
                  FREQ = rep( 0.5, 2 ),
                  SE   = rep( 0.00088, 2 ),
                  N    = rep( 485000, 2 ) )
k <- 0.008
df2 <- LmToOddsRatio( lmm.fil=df, k=k, std.err=1 )
df2$LOG_OR <- log(df2$OR_SE)
print(df2)

You can see that abs(logOR) is quite different even though the only thing that differs between the two SNPs is the sign of their LMM effect. Is there a bug in the code somewhere?

[lukelloydjones]

Hi kheilbron (sorry for not knowing real name),

I think this comes from the std.err=1 command, which returns the standard error of the OR not the OR itself. Set that flag to 0 and you will get back

 BETA FREQ      SE      N        OR

1 -0.0075 0.5 0.00088 485000 0.3589778 2 0.0075 0.5 0.00088 485000 2.7856877

which looks like what you want. If you want the SE of the OR we now recommend to do the process in the attached script.

Hope it helps,

Luke

From: kheilbron @.> Reply to: lukelloydjones/ORShiny @.> Date: Tuesday, 11 October 2022 at 12:25 pm To: lukelloydjones/ORShiny @.> Cc: Subscribed @.> Subject: [lukelloydjones/ORShiny] Discrepancy between negative and positive LMM effect conversions (Issue #1)

Hi Luke, thanks for an excellent paper and repo! I wanted to flag some odd behaviour I've been seeing with the function for converting LMM effects to odds ratios. Flipping the sign of the LMM effect should flip the sign of the resulting logOR, but the absolute value should be unchanged. However, I've found that positive and negative LMM effects get converted differently. Here is a minimal example for you:

df <- data.frame( BETA = c( -0.0075, 0.0075 ),

              FREQ = rep( 0.5, 2 ),

              SE   = rep( 0.00088, 2 ),

              N    = rep( 485000, 2 ) )

k <- 0.008

df2 <- LmToOddsRatio( lmm.fil=df, k=k, std.err=1 )

df2$LOG_OR <- log(df2$OR_SE)

print(df2)

You can see that abs(logOR) is quite different even though the only thing that differs between the two SNPs is the sign of their LMM effect. Is there a bug in the code somewhere?

— Reply to this email directly, view it on GitHubhttps://github.com/lukelloydjones/ORShiny/issues/1, or unsubscribehttps://github.com/notifications/unsubscribe-auth/ACL2HYZKUNSHDFSSWHG2PLTWCTF3FANCNFSM6AAAAAARB2WB2M. You are receiving this because you are subscribed to this thread.Message ID: @.***>

[karlheilbron]

Hi Luke,

Thank you so much for the speedy reply! You mentioned an attached script, but I'm not seeing one. Perhaps it didn't attach successfully?

In the meantime, my use case involves a GWAS Catalog dataset that does not provide FREQ, so I was hoping to use your script to estimate FREQ from SE, N, and BETA. If I'm not mistaken, lines 169-181 of your LmToOddsRatio() function calculate p0 and p1, which (in combination with k) be used to calculate FREQ. If I run just those lines using the BETA, SE, and N values from my minimal example:

beta <- c( -0.0075, 0.0075 )
se <- rep( 0.00088, 2 )
n  <- rep( 485000, 2 )
d  <- (beta ) / (2 * ((se ^ 2) * (n - 2) + (beta ^ 2)))
a  <- -beta
b  <- beta - 2 * d * k * beta
c  <- beta * 2 * k * (1 - k) * (d ^ 2) - k * (1 - k) * d + 
          k * beta * (d - (d^2))
p0_est.1 <- (-b - sqrt(b ^ 2 - 4 * a * c)) / (2 * a)
p0_est.2 <- (-b + sqrt(b ^ 2 - 4 * a * c)) / (2 * a)  
p0.prop  <- array(NA, length(p0_est.1))
p0.prop[which(beta < 0)]   <- p0_est.1[which(beta < 0)]
p0.prop[which(beta >= 0)]  <- p0_est.2[which(beta >= 0)]
p1_est <- p0.prop + d

And then if I calculate FREQ using:

FREQ = (1-k) * p0.prop + k * p1_est

I get FREQ = 0.0107 for both SNPs. However, FREQ_SNP1 should be equal to 1 - FREQ_SNP2 since we're modeling an allele flip here. This is making me think that there is an issue with calculating p0 and p1 from BETA, SE, and N. But perhaps you've already solved this in updated script that you mentioned?

Karl

[lukelloydjones]

Hi Karl,

I reviewed my last email and it had the attachment. It is small. Perhaps it’s because we are corresponding over the github service, which may block attachments but not sure.

I don’t think we can use the code to estimate the allele frequency. Been a long time since I worked on this and in this space. I think if you don’t have an allele frequency then you need to take one from an adequate reference population that aligns with your ancestry profile.

If you would like the SE calculation script then send me through a contact email and I will send again.

Hope this helps.

Luke

From: karlheilbron @.> Reply to: lukelloydjones/ORShiny @.> Date: Friday, 14 October 2022 at 8:23 am To: lukelloydjones/ORShiny @.> Cc: Mr Luke Lloyd-Jones @.>, Comment @.***> Subject: Re: [lukelloydjones/ORShiny] Discrepancy between negative and positive LMM effect conversions (Issue #1)

Hi Luke,

Thank you so much for the speedy reply! You mentioned an attached script, but I'm not seeing one. Perhaps it didn't attach successfully?

In the meantime, my use case involves a GWAS Catalog dataset that does not provide FREQ, so I was hoping to use your script to estimate FREQ from SE, N, and BETA. If I'm not mistaken, lines 169-181 of your LmToOddsRatio() function calculate p0 and p1, which (in combination with k) be used to calculate FREQ. If I run just those lines using the BETA, SE, and N values from my minimal example:

beta <- c( -0.0075, 0.0075 )

se <- rep( 0.00088, 2 )

n <- rep( 485000, 2 )

d <- (beta ) / (2 ((se ^ 2) (n - 2) + (beta ^ 2)))

a <- -beta

b <- beta - 2 d k * beta

c <- beta 2 k (1 - k) (d ^ 2) - k (1 - k) d +

      k * beta * (d - (d^2))

p0_est.1 <- (-b - sqrt(b ^ 2 - 4 a c)) / (2 * a)

p0_est.2 <- (-b + sqrt(b ^ 2 - 4 a c)) / (2 * a)

p0.prop <- array(NA, length(p0_est.1))

p0.prop[which(beta < 0)] <- p0_est.1[which(beta < 0)]

p0.prop[which(beta >= 0)] <- p0_est.2[which(beta >= 0)]

p1_est <- p0.prop + d

And then if I calculate FREQ using:

FREQ = (1-k) p0.prop + k p1_est

I get FREQ = 0.0107 for both SNPs. However, FREQ_SNP1 should be equal to 1 - FREQ_SNP2 since we're modeling an allele flip here. This is making me think that there is an issue with calculating p0 and p1 from BETA, SE, and N. But perhaps you've already solved this in updated script that you mentioned?

Karl

— Reply to this email directly, view it on GitHubhttps://github.com/lukelloydjones/ORShiny/issues/1#issuecomment-1278241039, or unsubscribehttps://github.com/notifications/unsubscribe-auth/ACL2HY5TMPLYG6GQSHLBULTWDCDU3ANCNFSM6AAAAAARB2WB2M. You are receiving this because you commented.Message ID: @.***>

[tom-a-bond]

Hi Luke, Would it be possible to add the script that you mentioned to calculate the SEs to this repos, or failing that to send me a copy please: tom [dot] bond [at] bristol [dot] ac [dot] uk Many thanks in advance Best wishes, Tom