How to model survival analysis when proportional hazards assumption is not met and stratification and time-varying are not possible? I am modelling a survival analysis over a rather long follow-up period (10 years). My exposure is time-invariant and clearly violates the proportional hazards assumptions so Cox Proportional Hazards regression models are not an option. I was wondering about alternatives to conduct my analyses. Please find below some key points:

*

*Stratification is not possible because it is my main variable of interest that violates the assumption and I need to compare between groups


*Time-varying models are not possible given the nature of my main variable of interest


*I initially thought about time-partitioned model (splitting follow-up time and interacting time with my main variable of interest) but I am not sure that is a good idea because when I plot the KM curves it is all crossing - so I struggle to find a good time interval for splitting
Given these considerations I have thought about employing flexible parametric models. Howevers, from my understanding they make strong assumptions about the shape of the curve - which is something I cannot be certain of. Would a flexible parametric model with restricted cubic spline what I am looking for? But How can I define the number of knots? And how about the distribution?
Could you please provide some inputs and examples? What do you suggest?
I use Stata MP 15 as statistical software
 A: I disagree that AFT and PO are necessarily the right next steps. It depends on what you are interested in learning from the model. If you are interested in estimating a hazard ratio, understand that, the idea that there is one hazard ratio that applies over a period of time, implies to some extent that hazards must be proportional.
On the other hand, in many applications there are much more informative summaries of survival analyses than hazard ratios, which don't inherently have a PH assumption baked in. For example, you can calculate risk differences and risk ratios at domain-relevant time-points. These are typically more intuitive and easier to interpret correctly than hazard ratios. RDs and RRs are still available using stratified Cox models (assuming your exposure is categorical). For an overview of these ideas, you can have a look at this reference:
https://pubmed.ncbi.nlm.nih.gov/25660080/
Now, if you insist on summarizing your data using hazard ratios, and hazards are not proportional, you can examine how the hazard ratio is changing over time using interactions between time and your time-invariant covariate of interest. This is a valid use of Cox models under non-proportional hazards and can be quite informative - explained in this paper: https://pubmed.ncbi.nlm.nih.gov/12915864/.  On the other hand, if the violation of proportionality is not too extreme, a single hazard ratio can still be a reasonable summary of the data - explain in this paper: https://pubmed.ncbi.nlm.nih.gov/32167523/.
A: Before you give up on proportional hazards (PH) modeling, make sure that you aren't relying on group-level Kaplan-Meier survival curves that haven't been adjusted fully for covariates. In particular, I wonder whether there might be some interactions with respect to survival between your main categorical variable of interest and other predictors that could be leading to these results. You seem to have enough data to do some fairly intricate modeling to evaluate those possibilities.
If PH still don't hold, then proportional odds (PO) or accelerated-failure-time (AFT) models would be next steps, as RaphaelS says in a comment. I don't know what's available in Stata for such modeling; otherwise, consider moving to R. The CRAN Survival Task View provides links to several implementations of such modeling approaches, which typically provide tools for examining goodness of fit appropriate to each model type.
Nevertheless, PO and AFT models have their own underlying assumptions. It's possible that your data won't adequately meet any of the PH, PO, or AFT assumptions. In that case, you might have to resort to stratification. That wouldn't necessarily be bad, it would just require more care in presenting results. For example, say that PH held for the other predictors in a Cox model stratified on your 6-level categorical exposure. Yes, there wouldn't be an association of the categorical exposures with event risk that held constant for all times. You could still, however, evaluate the associations of those other predictors with outcome, and compare outcomes among the exposure strata at specific times of interest.
