Survival models with competing risks - what alternatives to the cause-specific Cox models and Fine-Gray proportionate subdistribution hazards? I am working on a project involving survival analysis with various death causes. Specifically, there are 3 death causes and the question is if people with a particular medical history are more likely to die from one cause than another. There are other characteristics to be taken into account such as age, gender, etc. There seems to be two main approaches for the competing risk, cause-specific and Fine-Gray (https://onlinelibrary.wiley.com/doi/10.1002/sim.2712).
However, both approaches have some pitfalls - would anyone know good alternatives?
Machine-learning type, any, but those already realised in R or python function.
A bit more background: The cause-specific competing risk analysis is just a Cox model for each cause while others are deemed to be censoring events, and the hazard ratios describe the instantaneous risk to fail from that cause among those alive at that time. The pitfall is that cumulative incidence may not depend on a covariate in the same way, i.e. hazard ratio > 1 may not imply that expected share of people failed from this cause will be higher with the increase of that risk factor - as competing risk factors may skew incidence either way.
The Fine-Gray model directly imposes proportionate hazards on cumulative incidence functions, so it directly shows how a risk factor affects the incidence from that cause. However, it is not internally consistent and may be prone to bias as was published very recently in 2021: if one calculates event probabilities by each of the causes, the sum may be more than 1 (https://onlinelibrary.wiley.com/doi/full/10.1002/sim.9023)

I don't think you can circumnavigate these problems with alternative methodological solutions, since they stem from fundamental limitations in survival analysis.
Normally, in a Cox competing risks model death causes are considered independent from each other (i.e. censoring them when competing) because their degree of dependence is not directly measurable: you can observe only 1 death per patient. Nevertheless, working with cause-specific hazards allows you to obtain very punctual estimates of the impact different covariates have on your cause(s) of interest.
This is not the case when considering the cumulative incidence (which gives you the probability of occurrence of that cause, another quantity of interest). A covariate can have different effects on the hazard and cumulative incidence of the same cause. Intuitively, the reason is that the cumulative incidence of a given cause also depends on the cause-specific hazards of the competing causes. In your setting with 3 competing events, the cumulative incidence of cause 1 can be defined as$\ F_1 (t)= ∫_0^t[S(x)h_1 (x)dx] $, where$\ h_1 $ is the causes-specific hazard for cause 1 and $\ S $  is the overall survivor function. This is basically the probability of not dying due to any cause up to time$\ x$:$\ S(x)=exp \left\{ -∑_{i=1}^3 ∫_0^x[h_i (u)du] \right\} $.
The solution Fine & Gray proposed is to use a 'unnatural risk-set’, in which subjects exiting due to competing causes are kept as ‘at risk’ in the sample. This solves the issue above, but makes results rather tricky to interpret (or introduces other limitations, as the one you reported).
For a way out of this impasse, I would refer to this paper from Latouche, Allignol, Beyersmann, Labopin, and Fine: A competing risks analysis should report results on all cause-specific hazards and cumulative incidence functions, Journal of Clinical Epidemiology, (2013). The authors simply suggest to use all these approaches 'side-by-side', as pieces of the same puzzle, so that they complement each-other, allowing to build a complete understanding of the relationship of interest.
