Портал аспирантов - Параметрический анализ

Раз уж разбираем задачу, доводим её до логического абсурда. Предположим, что у нас распределение описывается некой вероятностной функцией, являющейся суммой двух вероятностных функций с весами, сумма которых равна единице (условие нормировки). Если вспомнить про интегралы, что интеграл суммы равен сумме интегралов, а постоянный член выносится за интеграл, то функции можно описать, например, так.

Код:

> pMyDist<-function(q,lambda1,lambda2,w=0.5) w*ppois(q,lambda1)+(1-w)*ppois(q,lambda2)

> dMyDist<-function(x,lambda1,lambda2,w=0.5) w*dpois(x,lambda1)+(1-w)*dpois(x,lambda2)

> qMyDist<-function(p,lambda1,lambda2,w=0.5) w*qpois(p,lambda1)+(1-w)*qpois(p,lambda2)

Где w принимает значение [0,1] и является весовым коэффициентом.
Если имеем p(k)=f(k,lambda), где f распределение Пуассона, то результирующая функция описывается как p(k)=w*f(k,lambda1)+(1-w)*f(k,lambda2)

Графически результат выглядит вот так (в сравнении с обычным распределением Пуассона) для функций
p₁(k)=f(k,7,2)
p₂(k)=0.8*f(k,7,2)+0.2*f(k,0,8)

Код:

> plot(1:25,dpois(1:25,7.2),xlab="x",ylab="f(x)")

> lines(dMyDist(1:25,7.2,4,0.8),type="p",col="red",pch=2)

> legend(x=15,y=0.1,c("Распределение Пуассона","Сложное распределение"),pch=c(1,2),col=c("black","red"))

http://aspirantura.spb.ru/forum/pict...&pictureid=986

Счастье заключается в том, что теперь эту функцию можно попробовать подогнать под наше распределение.

Код:

> XX<-fitdist(LT,"MyDist",start=list(lambda1=9 ,lambda2=5 ), method="mle")

> summary(XX)

Fitting of the distribution ' MyDist ' by maximum likelihood 

Parameters : 

        estimate Std. Error

lambda1 7.197402 0.05980069

lambda2 7.202613 0.05980405

Loglikelihood:  -123149.4   AIC:  246302.8   BIC:  246320.5 

Correlation matrix:

           lambda1    lambda2

lambda1  1.0000000 -0.9227544

lambda2 -0.9227544  1.0000000

Как видно, принципиально лучше не стало. Отрицательный результат -- тоже результат. Но метод может пригодиться. Я его как-то удачно использовал.