I'm trying to simulate a survival dataset, with the censoring variable $X$ following a continuous distribution with mass points. The distribution is $$f(x) = \begin{cases} \lambda_1e^{-\lambda_1x}, & 0 \le x \leq a \\ \lambda_2e^{-\lambda_2x}, & a < x \le b \\ \ \text{const}. & x > b \end{cases}$$ How to sample from such a distribution with $R$?
-
$\begingroup$ What is the distribution of $x$ $\endgroup$ – Dale C Nov 7 '20 at 7:47
-
$\begingroup$ @DaleC Sorry I forgot to mention that $X$ is the censoring variable. $\endgroup$ – Yujian Nov 7 '20 at 14:43
-
$\begingroup$ Right, but in order to sample from $f(x)$ you need to define the actual distribution of $x$ since your $f(x)$ is a compound distribution conditional on $x$. As in what is it's specific mass/density function? Otherwise how would you sample it? $\endgroup$ – Dale C Nov 7 '20 at 18:19
-
$\begingroup$ The variable X has a maximum right? Otherwise this pdf does not integrate to 1. $\endgroup$ – Sextus Empiricus Nov 7 '20 at 20:49
For your case it seems simplest to use inverse transform sampling.
Then you need to express the quantile function. For your case this requires I integrating your pdf which will be a piecewice exponential function and linear function. Then inverting this, which will be some function with logarithms.
$$F(x) = \begin{cases} 0 & x<0\\ c_1-e^{-\lambda_1x}, & 0 \le x \leq a \\ c_2-e^{-\lambda_2x}, & a < x \le b \\ c_3+x*const& b< x \le c_4 \\ 1 & x> c_4 \end{cases}$$
You will need to figure out those constants $c_i$ by setting the cases equal at the boundaries. E.g. for the last case $c_3 + x*const =1$ if $x=c_4$.
-
$\begingroup$ Thank you for pointing that out. Yes, I think it would be better to use CDF instead of PDF. I was thinking about using $b$ as a truncation point, so for any $x \geq b$, $F(x) = 1$. I should have expressed the $const$ as $1 - F(b)$, sorry about the confusion. Previously I was considering the inverse transform sampling (but I didn't find that page with the explanation). I thought it would only work for a continuous distribution, the CDF of which has a uniform distribution. Now it seems that the method also works for the distributions with mass points. That's great. Thank you. $\endgroup$ – Yujian Nov 7 '20 at 22:11
Just treat it as a conditional distribution. So sample $x$ according to its distribution, then nest a for loop for the appropriate bounds and append to a df.
df <- vector()
x = rdist(coefficients = [coef], size = 1)
n = int(size of trials)
count = 0
while (count < n){
if ( x > 0 & x<= a){
df = append(df, rpois(1,lamda1*x)
}
### repeat for the other bounds on x ###
count++
)
print(df)
-
$\begingroup$ Is 'rdist' the function for distance computation? $\endgroup$ – Yujian Nov 7 '20 at 14:50
-
1$\begingroup$ uhh sorry I was just generalising the generator function for the distribution of $x$ according to R notation for other distributions since you didn't state what distribution $x$ takes (rbinom etc.). It's technically an abuse of notation, but I mean, you have to have a defined distribution for $x$ if you're going to sample this compound distribution... $\endgroup$ – Dale C Nov 7 '20 at 18:16