0
$\begingroup$

If you run the following code, you will have a data frame real.dat which has 1063 samples for 20531 genes. There are 2 extra columns named time and event where time is the survival time and event is death in case of 1 and 0 in case of 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'

I want to get the univariate Cox regression p-value for each gene in the above data frame. Now, when I try to run the coxph function from the survival package even for one gene, it shows the following error -

> 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

What am I doing wrong here?

You can download the data from here.

$\endgroup$
1
$\begingroup$

I think the problem is that, for some reason, the genes in your data are not coded as numeric variables. For the HIF3A, for example:

> str(real.dat$HIF3A)
 Factor w/ 1059 levels "0.0000","0.3468",..: 309 327 294 108 611 458 347 346 576 417 ...
 - attr(*, "names")= chr [1:1063] "V2" "V3" "V4" "V5" ...

It should be numeric for the Cox model. The error you get is because coxph runs into trouble with fitting a model with a factor with 1059 levels. So yeah, check your data before fitting models.

NB: when converting factors to numeric, you should use as.numeric(as.character( )), as you can see from this:

> real.dat$HIF3A[1:5]
        V2     V3     V4     V5     V6 
    4.6261 4.7796 4.5313 2.2999 6.6477 
    1059 Levels: 0.0000 0.3468 0.4397 0.5484 0.5720 0.6396 0.9456 0.9635 0.9674 1.0976 1.1084 ... 9.9993
    > as.numeric(real.dat$HIF3A)[1:5]
[1] 309 327 294 108 611
> as.numeric(as.character(real.dat$HIF3A))[1:5]
[1] 4.6261 4.7796 4.5313 2.2999 6.6477
| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ Hi Theodor, many many thanks for taking the time to look at the data. That was indeed the problem here. Also, many thanks for suggesting the solution :). $\endgroup$ – nafizh Apr 20 '16 at 6:10

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.