jianyangqt / gcta

GCTA software
GNU General Public License v3.0
84 stars 26 forks source link

GCTA grm output not positive semi-definite #35

Closed willright28 closed 1 year ago

willright28 commented 1 year ago

Hi all,

I used the GCTA --make-grm function to get the binary and make the GRM using the following code:

a <- ReadGRMBin("test") m <- matrix(NA, ncol = length(a$diag), nrow = length(a$diag)) m[lower.tri(m)] <- a$off m[upper.tri(m)] <- t(m)[upper.tri(t(m))] diag(m) <- a$diag

When I used the R package rrBLUP function 'kin.blup', which requires a GRM input, received an error message "not positive semi-definite"

I've tried to add some noise to the diag but did not work. Is it relevant to my SNP dataset? For the GCTA input, I set missing rate, maf and LD threshold.

So could you give me any clue to solve this question? Thanks in advance!

jianyangqt commented 1 year ago

The GRM is non-positive definite because of losing one degree of freedom due to the use of the same genotype data for estimating allele frequencies. The GREML process is based on the phenotype variance-covariance matrix, which is often positive definite.

Best regards, Jian

Professor Jian Yang School of Life Sciences, Westlake University, Hangzhou Zhejiang 310024 China. T: +86-571-87089837. W: Lab websitehttps://yanglab.westlake.edu.cn/

From: willright28 @.> Reply-To: jianyangqt/gcta @.> Date: Friday, November 18, 2022 at 10:38 To: jianyangqt/gcta @.> Cc: Subscribed @.> Subject: [jianyangqt/gcta] GCTA grm output not positive semi-definite (Issue #35)

Hi all,

I used the GCTA --make-grm function to get the binary and make the GRM using the following code:

a <- ReadGRMBin("test") m <- matrix(NA, ncol = length(a$diag), nrow = length(a$diag)) m[lower.tri(m)] <- a$off m[upper.tri(m)] <- t(m)[upper.tri(t(m))] diag(m) <- a$diag

When I used the R package rrBLUP function 'kin.blup', which requires a GRM input, received an error message "not positive semi-definite"

I've tried to add some noise to the diag but did not work. Is it relevant to my SNP dataset? For the GCTA input, I set missing rate, maf and LD threshold.

So could you give me any clue to solve this question? Thanks in advance!

— Reply to this email directly, view it on GitHubhttps://apc01.safelinks.protection.outlook.com/?url=https%3A%2F%2Fgithub.com%2Fjianyangqt%2Fgcta%2Fissues%2F35&data=05%7C01%7Cjian.yang%40westlake.edu.cn%7C6b049f04dce24e11df7608dac90e07ef%7C7e82de2f7ef644169b9644c1457be81b%7C1%7C0%7C638043359344443858%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C3000%7C%7C%7C&sdata=BOpDQeHP3XDE%2BapOFOJAy9jAerhj%2BkypcuJdk6IQjAA%3D&reserved=0, or unsubscribehttps://apc01.safelinks.protection.outlook.com/?url=https%3A%2F%2Fgithub.com%2Fnotifications%2Funsubscribe-auth%2FAA4FPS3WSCY6APJSZZGFWPLWI3T3VANCNFSM6AAAAAASEAT6OY&data=05%7C01%7Cjian.yang%40westlake.edu.cn%7C6b049f04dce24e11df7608dac90e07ef%7C7e82de2f7ef644169b9644c1457be81b%7C1%7C0%7C638043359344443858%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C3000%7C%7C%7C&sdata=XmxahrSHF9TNKKnmaw%2FqQ973ckrvrHfkqn909Bk9neM%3D&reserved=0. You are receiving this because you are subscribed to this thread.Message ID: @.***>