Как выполнить анализ выживания для клинических данных?

Question

Мне нужно выполнить анализ выживаемости, чтобы найти значимые связи между конкретными путями деятельности и выживанием пациентов. Я пытаюсь выполнить анализ с помощью этого урока https://github.com/mforde84/RNAseq-Survival-Analysis-TCGA-KIRC/blob/master/survival_rnaseq_analysis.R. У меня есть два файла: клинические данные и uni_vals. Что-то не так, я думаю, в части анализа выживания. Сообщение выглядит следующим образом: «Ошибка во времени [[i]] <- сортировка (уникальная (у [кто, 1])): попытка выделить менее одного элемента в integerOneIndex "</p>

1) клинические данные

structure(list(Tumor_Sample_Barcode = structure(c(126L, 128L, 133L), .Label = c("TCGA-A7-A0CE", "TCGA-A7-A0CH", 
"TCGA-A7-A0DC"), class = "factor"), 
       classification_of_tumor = structure(c(1L, 1L, 1L), .Label   = "not reported", class = "factor"), 
last_known_disease_status = structure(c(1L, 1L,1L
), .Label = "not reported", class = "factor"), updated_datetime = structure(c(1L, 
1L, 1L), .Label = "2018-01-19T13:39:21.801433-06:00", class = "factor"), 
primary_diagnosis = structure(c(6L, 6L, 6L), .Label = c("C50.2", 
"C50.3", "C50.4", "C50.5", "C50.8", "C50.9", "C50.919"), class = "factor")), row.names = c(NA,5L), class = c("data.table", "data.frame"))

2) uni_vals

    structure(list(`TCGA-A7-A0CE` = c(0.800945270510658, 0.99887793401069, 
0.667341683672123, 0.999999999314536, 0.999999999314536), `TCGA-A7-A0CE.1` = c(0.778700980142054, 
0.998594728888895, 0.762898025094707, 0.999999999620033, 0.999999999620033
), `TCGA-A7-A0CH` = c(0.608118239725987, 0.992929539569249, 0.706256002082062, 
0.999999998256691, 0.999999998256691), `TCGA-A7-A0CH.1` = c(0.899309224249365, 
0.999869380713655, 0.797778413011216, 0.9999999997944, 0.9999999997944
), `TCGA-A7-A0DC` = c(0.535342987646728, 0.993464915142776, 0.409818699577936, 
0.999999996633627, 0.999999996633627)), row.names = c("Lipid degradation", 
"Lipid metabolism", "Chemotaxis", "Transcription regulation", 
"Transcription"), class = "data.frame")

3) мой код

получить индекс нормальных / контрольных образцов

  n_index <- which(substr(colnames(uni_vals),14,14) == '1')
    t_index <- which(substr(colnames(uni_vals),14,14) == '0')

ввод клинической информации

all_clin < -data.frame(cbind(clinical[,7],clinical[,10],clinical[,22]))
colnames(all_clin) <- c("new_tumor_days", "death_days", "followUp_days")
rownames(all_clin) <- clinical$Tumor_Sample_Barcode

время до опухоли

all_clin$new_time <- c() for (i in 1:length(as.numeric(as.character(all_clin$new_tumor_days)))){all_clin$new_time[i] <-ifelse(is.na(as.numeric(as.character(all_clin$new_tumor_days))[i]), 
                            as.numeric(as.character(all_clin$followUp_days))[i],
                            as.numeric(as.character(all_clin$new_tumor_days))[i])}

время до смерти

all_clin$new_death <- c()
for (i in 1:length(as.numeric(as.character(all_clin$death_days)))){
  all_clin$new_death[i] <- ifelse(is.na(as.numeric(as.character(all_clin$death_days))[i]),
                                  as.numeric(as.character(all_clin$followUp_days))[i],
                                  as.numeric(as.character(all_clin$death_days))[i])
}

событие цензора смерти

 all_clin$death_event <- ifelse(clinical$vital == "alive", 0, 1)

фильтрации

colnames(uni_vals) <- gsub("\\.","-",substr(colnames(uni_vals),1,12))

соответствует all_clin и uni_vals

 ind_tum <- which(unique(colnames(uni_vals)) %in% rownames(all_clin))
    ind_clin <- which(rownames(all_clin) %in% colnames(uni_vals))

функция интереса

ind_func <- which(rownames(uni_vals) == "Lipid degradation")

создать вектор событий для унифицированных данных

event_uni <- t(apply(uni_vals, 1, function(x) ifelse(abs(x) <0.01,1,0)))

анализ выживания

  s <- survfit(Surv(as.numeric(as.character(all_clin$new_death))[ind_clin], 
                  all_clin$death_event[ind_clin]) ~ event_uni[ind_func, ind_tum])
    s1 <- tryCatch( survdiff(Surv(as.numeric(as.character(all_clin$new_death))[ind_clin], 
                all_clin$death_event[ind_clin]) ~ event_uni[ind_func, ind_tum]), 
      error = function(e) 
    return(NA))

Answer 1

Я не уверен, в чем причина вашей проблемы, но ваш код выдает ошибку уже в первой части «1) Клинические данные»: «Ошибка в as.character.factor (x): неправильный коэффициент».Я бы посоветовал сначала проверить это - если эта часть и структуры данных неверны, последующие блоки кода и анализ выживания не будут работать.

Также у вас есть опечатка в:

##input clinical information
all_clin < -data.frame(cbind(clinical[,7],clinical[,10],clinical[,22]))

Как должно быть:

##input clinical information
all_clin <- data.frame(cbind(clinical[,7],clinical[,10],clinical[,22]))

Если это не помогло, но ошибка не устранена, отметьте этот пост , где описано решение для подобной ошибки, и дайте мне знать.