Как получить p-значение Cox для каждого гена?

Question

Если вы запустите следующий код, у вас будет кадр данных real.dat, который имеет 1063 образца для 20531 генов. Есть 2 дополнительные столбцы с именем time и event где time это время выживания и event является death в случае 1 и 0 в случае censored.

lung.dat <- read.table("genomicMatrix_lung") 
lung.clin.dat <- read.delim("clinical_data_lung") 

# For clinical data, get only rows which do not have NA in column "X_EVENT" 
lung.no.na.dat <- lung.clin.dat[!is.na(lung.clin.dat$X_EVENT), ] 

# Getting the transpose of main lung cancer data 
ge <- t(lung.dat) 
# Getting a vector of all the id's in the clinical data frame without any 'NA' values 
keep <- lung.no.na.dat$sampleID 

# getting only the samples(persons) for which we have a value rather than 'NA' values 
real.dat <- ge[ge[, 1] %in% keep, ] 

# adding the 2 columns from clinical data to gene expression data 
keep_again <- real.dat[, 1] 
temp_df <- lung.no.na.dat[lung.no.na.dat$sampleID %in% keep_again, ] 

# naming the columns into our gene expression data 
col_names <- ge[1, ] 
colnames(real.dat) <- col_names 

dd <- temp_df[, c('X_TIME_TO_EVENT', 'X_EVENT')] 
real.dat <- cbind(real.dat, dd) 

# renaming the 2 new added columns 
colnames(real.dat)[colnames(real.dat) == 'X_TIME_TO_EVENT'] <- 'time' 
colnames(real.dat)[colnames(real.dat) == 'X_EVENT'] <- 'event' 

Я хочу получить одномерную регрессионную p-значение Cox для каждого гена в вышеуказанном кадре данных. Как я могу это получить?

Вы можете скачать данные с here.

Edit: Извините, что недостаточно разъясняю. Я уже пытался получить его с помощью функции coxph из библиотеки survival. Но даже для одного гена, он показывает следующее сообщение об ошибке -

> coxph(Surv(time, event) ~ HIF3A, real.dat) Error in fitter(X, Y, strats, offset, init, control, weights = weights, : NA/NaN/Inf in foreign function call (arg 6) In addition: Warning message: In fitter(X, Y, strats, offset, init, control, weights = weights, : Ran out of iterations and did not converge

Вот почему я не обеспечивал меньший воспроизводимый пример.

Answer 1

Вы действительно будете делать одномерную регрессию для каждого гена 20531 генов?

Гадание дико на структуру данных (так, создавая фиктивный набор, основанный на примерах в помощи), и угадать, что вы пытаетесь сделать на следующий примере игрушки .....

library("survival") 
?coxph ## to see the examples 
## create dummy data 
test <- list(time=c(4,3,1,1,2,2,3), 
      event=c(1,1,1,0,1,1,0), 
      gene1=c(0,2,1,1,1,0,0), 
      gene2=c(0,0,0,0,1,1,1)) 
## Cox PH regression 
coxph(Surv(time, event) ~ gene1, test) 
coxph(Surv(time, event) ~ gene2, test) 

Возможно, вы захотите использовать следующее, чтобы получить CI и дополнительную информацию.

summary(coxph(...)) 

Надеется, что код достаточно воспроизводимый, чтобы помочь вам прояснить вопрос