-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathNPHpdf.R
51 lines (40 loc) · 1.68 KB
/
NPHpdf.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
"Esta función calcula la densidad de una distribución
NPH. Con las variables de salidad la densidad en si.
f densidad
Abajo se definen las variables del modelo
y es un vector de puntos de datos
alpha es le valor inicial de la distribución PH
T0 es la matrix subintensidad de la PH
theta contiene los parametros de la distribución de
escala.
SCANLIGCDF contiene el nombre de la distribución de
escala.
El soporte discreto sobre el que se discretizará la
distribución a escala."
nphpdf<-function(y,MODEL,PARAMETERS){
SCALINGCDF <- MODEL$SCALINGCDF
S <- MODEL$S
alpha <- PARAMETERS$alpha
T0 <- PARAMETERS$T0
theta <- PARAMETERS$theta
PI <- nphpmf(SCALINGCDF,S,theta)
PI$PI
"%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%% Technical stuff (size of vectors, initializing variables) %%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%"
Ny=length(y) # Ny = number of data points
NPi=length(PI ) # Npi = number of levels
f=matrix(rep(0,Ny),1) # Initializing vector of densities
t=apply(T0,1,sum)*(-1) # Absorbtion rates
"%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%% CALCULATION OF THE DENSITY %%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% The index i corresponds to the levels (rows of the matrix fy)
% while the index j corresponds to data points (columns)"
for(j in 1:Ny){
for(i in 1:NPi){
f[j]<-f[j]+PI[i]%*%alpha%*%(expm(T0*y[j]/S[i])+1)*t/S[i]
}
}
return(f)
}