Фиксирование кривых плотности выживания с использованием разных распределений

Question

Я работаю с некоторыми нормальными данными, и, естественно, хочу продемонстрировать результаты логарифмического распределения в лучшем перекрытии, чем другие возможные распределения. По существу, я хочу повторить следующий график с моими данными:

где кривыми подогнанные плотностей сопоставляются над log(time).

текста, где связанное изображение от описывает процесс установки каждую модели и получению следующих параметров:

Для этой цели, я установлен четыре наивных моделей выживания с вышеупомянутыми распределениями :

survreg(Surv(time,event)~1,dist="family") 

и экстрагируют параметр формы (α) и коэффициента (β).

У меня есть несколько вопросов, касающихся процесса:

1) Является ли это правильный путь идти об этом? Я просмотрел несколько R-пакетов, но не смог найти ту, что графиков плотности графиков, как встроенная функция, поэтому я чувствую, что я должен упускать из виду что-то очевидное.

2) Имеют ли значения соответствующие логарифмически нормальное распределение (μ и σ $^2 $) только среднее и дисперсия перехвата?

3) Как создать подобную таблицу в R? (Может быть, это скорее вопрос переполнения стека). Я знаю, что могу просто cbind их вручную, но меня больше интересует их вызов из встроенных моделей. survreg объекты хранят оценки коэффициентов, но при вызове survreg.obj$coefficients получается именованный вектор числа (а не просто число).

4) Самое главное, как я могу построить аналогичный граф? Я думал, что это будет довольно просто, если я просто извлечу параметры и нарисую их над гистограммой, но пока не повезло. Автор текста говорит, что он оценил кривые плотности по параметрам, но я просто получаю точечную оценку - чего мне не хватает? Должен ли я вычислять кривые плотности вручную на основе распределения перед построением графика?

Я не уверен, как обеспечить mwe в этом случае, но, честно говоря, мне просто нужно общее решение для добавления нескольких кривых плотности для данных о выживании. С другой стороны, если вы считаете, что это поможет, не стесняйтесь рекомендовать решение mwe, и я постараюсь его создать.

Спасибо за ваш ввод!

Редактировать: На основании сообщения eclark я сделал некоторый прогресс. Мои параметры:

Dist = data.frame(
Exponential = rweibull(n = 10000, shape = 1, scale = 6.636684), 
Weibull = rweibull(n = 10000, shape = 6.068786, scale = 2.002165), 
Gamma = rgamma(n = 10000, shape = 768.1476, scale = 1433.986), 
LogNormal = rlnorm(n = 10000, meanlog = 4.986, sdlog = .877) 
) 

Однако, учитывая огромную разницу в масштабах, это то, что я получаю:

Возвращаясь к вопросу № 3, это, как я должен получить параметры ? В настоящее время это, как я делаю это (извините за беспорядок):

summary(fit.exp) 

Call: 
survreg(formula = Surv(duration, confterm) ~ 1, data = data.na, 
dist = "exponential") 
     Value Std. Error z p 
(Intercept) 6.64  0.052 128 0 

Scale fixed at 1 

Exponential distribution 
Loglik(model)= -2825.6 Loglik(intercept only)= -2825.6 
Number of Newton-Raphson Iterations: 6 
n= 397 

summary(fit.wei) 

Call: 
survreg(formula = Surv(duration, confterm) ~ 1, data = data.na, 
dist = "weibull") 
     Value Std. Error z  p 
(Intercept) 6.069  0.1075 56.5 0.00e+00 
Log(scale) 0.694  0.0411 16.9 6.99e-64 

Scale= 2 

Weibull distribution 
Loglik(model)= -2622.2 Loglik(intercept only)= -2622.2 
Number of Newton-Raphson Iterations: 6 
n= 397 

summary(fit.gau) 

Call: 
survreg(formula = Surv(duration, confterm) ~ 1, data = data.na, 
dist = "gaussian") 
     Value Std. Error  z  p 
(Intercept) 768.15 72.6174 10.6 3.77e-26 
Log(scale) 7.27  0.0372 195.4 0.00e+00 

Scale= 1434 

Gaussian distribution 
Loglik(model)= -3243.7 Loglik(intercept only)= -3243.7 
Number of Newton-Raphson Iterations: 4 
n= 397 

summary(fit.log) 

Call: 
survreg(formula = Surv(duration, confterm) ~ 1, data = data.na, 
dist = "lognormal") 
     Value Std. Error z   p 
(Intercept) 4.986  0.1216 41.0 0.00e+00 
Log(scale) 0.877  0.0373 23.5 1.71e-122 

Scale= 2.4 

Log Normal distribution 
Loglik(model)= -2624 Loglik(intercept only)= -2624 
Number of Newton-Raphson Iterations: 5 
n= 397 

Я чувствую, что я особенно портя Логнормальное, учитывая, что это не стандартная форма-и-коэффициент тандем, но среднее и дисперсия.

Answer 1

Попробуйте это; идея генерации случайных переменных с помощью случайных функций distribtion, а затем построение функции плотности с выходными данными, вот пример, как вам нужно:

require(ggplot2) 
require(dplyr) 
require(tidyr) 

SampleData <- data.frame(Duration=rlnorm(n = 184,meanlog = 2.859,sdlog = .246)) #Asume this is data we have sampled from a lognormal distribution 

#Then we estimate the parameters for different types of distributions for that sample data and come up for this parameters 
#We then generate a dataframe with those distributions and parameters 
Dist = data.frame(
    Weibull = rweibull(10000,shape = 1.995,scale = 22.386), 
    Gamma = rgamma(n = 10000,shape = 4.203,scale = 4.699), 
    LogNormal = rlnorm(n = 10000,meanlog = 2.859,sdlog = .246) 
) 

#We use gather to prepare the distribution data in a manner better suited for group plotting in ggplot2 
Dist <- Dist %>% gather(Distribution,Duration) 

#Create the plot that sample data as a histogram 
G1 <- ggplot(SampleData,aes(x=Duration)) + geom_histogram(aes(,y=..density..),binwidth=5, colour="black", fill="white") 

#Add the density distributions of the different distributions with the estimated parameters 
G2 <- G1 + geom_density(aes(x=Duration,color=Distribution),data=Dist) 

plot(G2)