9
$\begingroup$

I am trying to run a Cox regression on a sample 2,000,000 row dataset as follows using only R. This is a direct translation of a PHREG in SAS. The sample is representative of the structure of the original dataset.

##
library(survival)

### Replace 100000 by 2,000,000

test <- data.frame(start=runif(100000,1,100), stop=runif(100000,101,300), censor=round(runif(100000,0,1)), testfactor=round(runif(100000,1,11)))

test$testfactorf <- as.factor(test$testfactor)
summ <- coxph(Surv(start,stop,censor) ~ relevel(testfactorf, 2), test)

# summary(summ)
##

user  system elapsed 
9.400   0.090   9.481

The main challenge is in the compute time for the original dataset (2m rows). As far as I understand, in SAS this could take up to 1 day, ... but at least it finishes.

  • Running the example with only 100,000 observations take only 9 seconds. Thereafter the time increases almost quadratically for every 100,000 increment in the number of observations.

  • I have not found any means to parallelize the operation (e.g., we can leverage a 48-core machine if this was possible)

  • Neither biglm nor any package from Revolution Analytics is available for Cox regression, and so I cannot leverage those.

Is there a means to represent this in terms of a logistic regression (for which there are packages in Revolution) or if there are any other alternatives to this problem? I know that they are fundamentally different, but it's the closest I can assume as a possibility given the circumstances.

Improve this question
$\endgroup$
2
  • 2
    $\begingroup$ Conditional logistic regression and Cox regression are strictly related. stayconsistent.wordpress.com/2009/01/25/… $\endgroup$
    – boscovich
    Commented Jun 16, 2013 at 8:36
  • $\begingroup$ were you able to arrive at an elegant solution to this? I am calculating HRs from $coxph$ for thousands of genes over multiple dataets, and $coxph$ is a bottle neck. At present I am using $apply$ to loop over genes. PS stack is not letting me comment as I am new user. $\endgroup$
    – Arshi Arora
    Commented Mar 12, 2019 at 21:03

2 Answers 2

Reset to default
3
$\begingroup$

I run cox regression on a 7'000'000 observation dataset using R and this is not a problem. Indeed, on bivariate models I get the estimates in 52 seconds. I suggest that it is -as often with R- a problem related to the RAM available. You may need at least 12GB to run the model smoothly.

Improve this answer
$\endgroup$
1
$\begingroup$

I went directly to the hardcore fit function (agreg.fit), which under the hood is called for the computations:

n <- nrow(test)
y <- as.matrix(test[, 1:3])
attr(y, "type") <- "right"
x <- matrix(1:11, n, 11, byrow=TRUE)
colnames(x) <- paste("level", 1:11, sep="")
x <- x[, -2] == test$testfactor
mode(x) = "numeric"

fit2 <- agreg.fit(x, y, strata=NULL, control=coxph.control(), method="efron", 
                  init=rep(0, 10), rownames=1:n)

However, the time elapsed when doubling the sample size gets quadratic as you mentioned. Also decreasing the epsilon in coxph.control does not help.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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