One of the central goals of precision health is the understanding and interpretation of high-dimensional biological data to identify genes and markers associated with disease initiation, development and outcomes. Significant effort has been committed to harness gene expression data as real-valued matrices for multiple analyses while accounting for time-to-event modeling by including survival times. Traditional biological analysis has focused separately on non-negative matrix factorization (NMF) of the gene expression data matrix and survival regression with Cox proportional hazards model. In this work, Cox proportional hazards regression is integrated with NMF by imposing survival constraints. This is accomplished by jointly optimizing the Frobenius norm and partial log likelihood for events such as death or relapse. Simulation results based on synthetic data demonstrated the superiority of the proposed methodology, when compared to other NMF algorithms, in finding survival associated gene clusters. In addition, using breast cancer gene expression data, the proposed technique can unravel critical clusters of cancer genes. The discovered gene clusters reflect rich biological implications and can help identify survival-related biomarkers. Towards the goal of precision health and cancer treatments, the proposed algorithm can help understand and interpret high-dimensional heterogeneous genomics data with accurate identification of survival-associated gene clusters.