2017-05-01 1 views
1

Ich möchte in der Lage sein, die p-Werte für die lmekin-Objekt von der Coxme-Paket produziert zu sehen.Wie p-Werte aus lmekin Objekte in Coxme-Paket extrahieren

z.

model= lmekin(formula = height ~ score + sex + age + (1 | IID), data = phenotype_df, varlist = kinship_matrix) 

Ich habe versucht:

summary(model) 
coef(summary(model)) 
summary(model$coefficient$fixed) 
fixef(model)/ sqrt(diag(vcov(model)) #(Calculates Z-scores but not p-values) 

Aber diese hat nicht funktioniert. Wie sehe ich die p-Werte für dieses lineare gemischte Modell?

Antwort

2

Es dauerte Jahre der Suche, um dies herauszufinden, aber ich bemerkte eine Menge anderer ähnlicher Fragen ohne richtige Antworten, also wollte ich das beantworten.

Sie verwenden:

library(coxme) 
print(model) 
  • Hinweis Es ist wichtig, das coxme Paket vorab zu laden oder es wird nicht funktionieren.

Ich habe bemerkt, auch viele Beiträge, wie der p-Wert von lmekin Objekten zu extrahieren, oder wie der p-Wert von coxme Objekten im Allgemeinen zu extrahieren. Ich habe diese Funktion geschrieben, die auf dem Funktionscode coxme:::print.coxme basiert (um den Codetyp coxme:::print.coxme direkt in R anzuzeigen). print berechnet p-Werte im laufenden Betrieb - diese Funktion ermöglicht die Extraktion von p-Werten und speichert sie in einem Objekt.

Beachten Sie, dass mod Ihr Modell ist, z. mod <- lmekin(y~x+a+b) Verwenden Sie print(mod), um zu überprüfen, ob die Tabellen übereinstimmen.

extract_coxme_table <- function (mod){ 
    beta <- mod$coefficients$fixed 
    nvar <- length(beta) 
    nfrail <- nrow(mod$var) - nvar 
    se <- sqrt(diag(mod$var)[nfrail + 1:nvar]) 
    z<- round(beta/se, 2) 
    p<- signif(1 - pchisq((beta/se)^2, 1), 2) 
    table=data.frame(cbind(beta,se,z,p)) 
    return(table) 
} 
+0

Ich musste die erste Zeile der Funktion durch 'beta <- fixef (mod)' ersetzen. – Axeman

Verwandte Themen