A template for the R scripts, that were used to generate the data sets and analyze the estimation results for full models, is given in the following. For the other two model types, the script was slightly changed. The values for the variables n, g and sigma where chosen depending on the sample size type and the correlation structure of the simulation. Further, the seed was varied across the simulations.
# L o a d p a c k a g e s l i b r a r y( s t r i n g r ) l i b r a r y( t i d y v e r s e ) l i b r a r y( g a m l s s ) l i b r a r y( l a t t i c e ) l i b r a r y( P e n C o x F r a i l ) l i b r a r y( m v t n o r m )
#
-# 1. S i m u l a t i o n Set - up
# For the s i m u l a t i o n , n i n d i v i d u a l s are a s s u m e d . For e a c h i n d i v i d u a l , n _ i o b s e r v a t i o n s are g e n e r a t e d f r o m the DGP , r e s u l t i n g in N o b s e r v a t i o n s .
#
-# N u m b e r s of o b s e r v a t i o n s and g r o u p s n < - n
g < - g
# I n f l u e n c e on the p a r a m e t e r s
f . mu < - f u n c t i o n( x ) {2 + 0.2 * x + sin( x ) /4}
f . s i g m a < - f u n c t i o n( x ) { -( x -5) ^2 / 20 - 1}
f . nu < - f u n c t i o n( x ) { - 0.5 + ( - 0 . 2 5 ) * x }
# L i n k F u n c t i o n i n c l . r a n d o m i n t e r c e p s
mu . l i n k < - f u n c t i o n( x , RI ) {exp( f . mu ( x ) + RI ) }
s i g m a . l i n k < - f u n c t i o n( x , RI ) {exp( f . s i g m a ( x ) + RI ) } nu . l i n k < - f u n c t i o n( x , RI ) { 1 / ( 1 +exp( - ( f . nu ( x ) + RI ) ) ) }
# G e n e r a t e R a n d o m I n t e r c e p t s s e t . s e e d( s e e d )
# C o v a r i a n c e m a t r i x s i g m a < - s i g m a
# D r a w r a n d o m i n t e r c e p t s
RI < - r m v n o r m ( g , m e a n = c(0 ,0 ,0) , s i g m a = s i g m a )
# C e n t e r and s c a l e to m e a n and v a r i a n c e RI < - t(t( RI ) - c o l M e a n s( RI ) )
RI [ ,1] < - RI [ ,1]/sd( RI [ ,1]) * s q r t( 0 . 1 ) RI [ ,2] < - RI [ ,2]/sd( RI [ ,2]) * s q r t( 0 . 5 ) RI [ ,3] < - RI [ ,3]/sd( RI [ ,3]) * s q r t( 0 . 2 5 ) r o u n d(c o l M e a n s( RI ) , 6)
# G e n e r a t e d a t a f r a m e w i t h x and r a n d o m i n t e r c e p t s df < - d a t a.f r a m e( x = r u n i f( n , min = 0 , max = 10) , g r o u p s = 1: g ,
R I m u = RI [ ,1] , R I s i g m a = RI [ ,2] , R I n u = RI [ ,3]) df $ g r o u p s < - as.f a c t o r( df $ g r o u p s )
# G e n e r a t e d a t a f r a m e w i t h c o r r e s p o n d i n g mu , s i g m a and nu mu < - mu . l i n k ( df $ x , df $ R I m u )
s i g m a < - s i g m a . l i n k ( df $ x , df $ R I s i g m a ) nu < - nu . l i n k ( df $ x , df $ R I n u )
df $ yn < - r Z A G A( n , mu = mu , s i g m a = sigma , nu = nu )
# H e l p i n g v a r i a b l e
df $ y < - i f e l s e( df $ yn == 0 , 0 , 1)
# S a v e the d a t a f r a m e to c a l c u l a t e the m o d e l in B a y e s X : d f B a y e s X < - df [ , c(" x ", " g r o u p s ", " yn ", " y ") ]
w r i t e.t a b l e( d f B a y e s X , p a s t e 0(g e t w d() , " / d a t a _ ",n ," _ ",g ," . raw ") , q u o t e = FALSE, row.n a m e s = F A L S E)
#
-# 2. V i s u a l i z a t i o n
# V i s u a l i z e the t r u e i n f l u e n c e of x on the d i s t r i b u t i o n a l
p a r a m e t e r s g r o u p e d by the r a n d o m i n t e r c e p t set to the 0.01 , 0.1 , 0.25 , 0.5 , 0.75 , 0.9 and 0 . 9 5 q u a n t i l e of the u n d e r l y i n g n o r m a l d i s t r i b u t i o n s
#
-# T h e o r e t i c a l r a n d o m i n t e r c e p t s
RI . t h e o < - m a t r i x(0 , n r o w = 7 , n c o l = 3)
RI . t h e o [ ,1] < - q n o r m(c(0.95 , 0 . 9 , 0 . 7 5 , 0 . 5 , 0 . 2 5 , 0 . 1 , 0 . 0 5 ) , m e a n = 0 , sd = s q r t( 0 . 1 ) )
RI . t h e o [ ,2] < - q n o r m(c(0.95 , 0 . 9 , 0 . 7 5 , 0 . 5 , 0 . 2 5 , 0 . 1 , 0 . 0 5 ) , m e a n = 0 , sd = s q r t( 0 . 5 ) )
RI . t h e o [ ,3] < - q n o r m(c(0.95 , 0 . 9 , 0 . 7 5 , 0 . 5 , 0 . 2 5 , 0 . 1 , 0 . 0 5 ) , m e a n = 0 , sd = s q r t( 0 . 2 5 ) )
# E f f e c t of x on mu xx < - seq( 0 . 1 , 1 0 , 0 . 0 1 )
q u a n t < - p a s t e(c(0.95 , 0 . 9 , 0 . 7 5 , 0 . 5 , 0 . 2 5 , 0 . 1 , 0 . 0 5 ) )
col < - c o l o r R a m p P a l e t t e (c(" d e e p s k y b l u e 2 ", " b l u e 2 ", " b l u e 4 ", "
b l a c k ", " b l u e 4 "," b l u e 2 ", " d e e p s k y b l u e 2 ") ) (7) lwd < - c(2 ,2 ,2 ,4 ,2 ,2 ,2)
par ( m f r o w = c(1 ,3) , mar = c(7.1 , 6.1 , 8.1 , 6 . 1 ) )
p l o t( NA , y l a b = e x p r e s s i o n(p a s t e( mu ," ( x ) ") ) , m a i n = " ", y l i m = c (0 ,100) , lwd = 3 , cex . lab = 3 , x a x t = ’ n ’, y a x t = ’ n ’, x l a b = "
x ", x l i m = c( 0 , 1 0 . 6 ) )
t i t l e(e x p r e s s i o n(a t o p(p a s t e(" E f f e c t of ", x , " on ", mu , "
g r o u p e d by ") , p a s t e(" d i f f e r e n t q u a n t i l e s for the r a n d o m i n t e r c e p t ") ) ) , cex . m a i n = 3)
a x i s(1 ,c e x . a x i s=2) a x i s(2 ,c e x . a x i s=2) j < - 1
a b s t a n d < - c( -0.5 , 0.5 ,0 ,0 ,0 ,0 ,0)
lty < - c(" d o t t e d ", " d o t t e d ", " d o t t e d ", " s o l i d ", " d a s h e d ", " d a s h e d
", " d a s h e d ")
for( i in RI . t h e o [ ,1]) {
l i n e s( xx , mu . l i n k ( xx , i ) , col = col [ j ] , lwd = lwd [ j ] , lty = lty [ j ])
t e x t(10.5 , mu . l i n k (10 , i ) + a b s t a n d [ j ] , l a b e l s = q u a n t [ j ] , col = col [ j ] , cex = 2)
j < - j + 1 }
l e g e n d(" t o p l e f t ", l e g e n d = c(s p r i n t f(" %.3 f ", RI . t h e o [ 1 : 3 , 1 ] ) , 0 , s p r i n t f(" %.3 f ", RI . t h e o [ 5 : 7 , 1 ] ) ) , col = col , lty = lty , lwd = c(3 ,3 ,3 ,4 ,3 ,3 ,3) , cex = 1.8 , t i t l e = " C o n c r e t e v a l u e s ")
# E f f e c t of x on s i g m a
p l o t( NA , y l a b = e x p r e s s i o n(p a s t e( sigma ," ( x ) ") ) , y l i m = c(0 ,1.5) , m a i n = " ", lwd = 3 , cex . lab = 3 , x a x t = ’ n ’, y a x t = ’ n ’, x l a b
= " x ", x l i m = c(0 ,10) )
t i t l e(e x p r e s s i o n(a t o p(p a s t e(" E f f e c t of ", x , " on ", sigma , "
g r o u p e d by ") , p a s t e(" d i f f e r e n t q u a n t i l e s for the r a n d o m i n t e r c e p t ") ) ) , cex . m a i n = 3)
a x i s(1 ,c e x . a x i s=2) a x i s(2 ,c e x . a x i s=2) j < - 1
a b s t a n d < - c( 0 . 0 4 , 0 . 0 4 , 0 . 0 4 , 0 . 0 4 , 0 . 0 4 , 0 . 0 4 , - 0 . 0 4 ) for( i in RI . t h e o [ ,2]) {
l i n e s( xx , s i g m a . l i n k ( xx , i ) , col = col [ j ] , lwd = lwd [ j ] , lty
= lty [ j ])
t e x t(5 , s i g m a . l i n k (5 , i ) + a b s t a n d [ j ] , l a b e l s = q u a n t [ j ] , col
= col [ j ] , cex = 2) j < - j + 1
}
l e g e n d(" t o p l e f t ", l e g e n d = c(s p r i n t f(" %.3 f ", RI . t h e o [ 1 : 3 , 2 ] ) , 0 , s p r i n t f(" %.3 f ", RI . t h e o [ 5 : 7 , 2 ] ) ) , col = col , lty = lty , lwd = c(3 ,3 ,3 ,4 ,3 ,3 ,3) , cex = 1.8 , t i t l e = " C o n c r e t e v a l u e s ")
# E f f e c t of x on nu
p l o t( NA , y l a b = e x p r e s s i o n(p a s t e( nu ," ( x ) ") ) , y l i m = c(0 ,0.6) , m a i n = " ", lwd = 3 , cex . lab = 3 , x a x t = ’ n ’, y a x t = ’ n ’, x l a b =
" x ", x l i m = c( -0.6 ,10) )
t i t l e(e x p r e s s i o n(a t o p(p a s t e(" E f f e c t of ", x , " on ", nu , "
g r o u p e d by ") , p a s t e(" d i f f e r e n t q u a n t i l e s for the r a n d o m i n t e r c e p t ") ) ) , cex . m a i n = 3)
a x i s(1 ,c e x . a x i s=2) a x i s(2 ,c e x . a x i s=2) j < - 1
a b s t a n d < - c(0.01 ,0.01 , 0 . 0 1 , 0 . 0 1 , 0 . 0 0 5 , 0 , - 0 . 0 3 ) for( i in RI . t h e o [ ,3]) {
l i n e s( xx , nu . l i n k ( xx , i ) , col = col [ j ] , lwd = lwd [ j ] , lty = lty [ j ])
t e x t( -0.5 , nu . l i n k (0 , i ) , l a b e l s = q u a n t [ j ] , col = col [ j ] , cex = 2)
j < - j + 1 }
l e g e n d(" t o p r i g h t ", l e g e n d = c(s p r i n t f(" %.3 f ", RI . t h e o [ 1 : 3 , 3 ] ) , 0 , s p r i n t f(" %.3 f ", RI . t h e o [ 5 : 7 , 3 ] ) ) , col = col , lty = lty , lwd = c(3 ,3 ,3 ,4 ,3 ,3 ,3) , cex = 1.8 , t i t l e = " C o n c r e t e v a l u e s ")
#
-# 3. R e s u l t s
# It is v e r y i m p o r t a n t to m e n t i o n t h a t the p a r a m e t r i z a t i o n t h a t was u s e d to g e n e r a t e the d a t a are d i f f e r e n t for s i g m a and nu to the p a r a m e t r i z a t i o n in B a y e s X :
#
# * mu : The s a m e p a r a m e t e r i z a t i o n is u s e d .
# * s i g m a : I h a v e g e n e r a t e d the i n f l u e n c e of x on the p a r a m e t e r s i g m a ( as you can see a b o v e ) . In c o n t r a s t , i n s t e a d of sigma , B a y e s X u s e s tau for the s e c o n d p a r a m e t e r . The c o n n e c t i o n b e t w e e n tau and s i g m a is tau = 1/ s i g m a ^2
# * nu : I h a v e m o d e l e d the i n f l u e n c e of x on the p a r a m e t e r nu . In BayesX , the p o i n t m a s s at z e r o is m o d e l e d by pi w h i c h is pi = 1 - nu
# In o r d e r to a c c o u n t for the d i f f e r e n t p a r a m e t r i z a t i o n , the e s t i m a t e s are t r a n s f o r m e d to o b t a i n the e s t i m a t e d e f f e c t s on s i g m a and nu to c o m p a r e t h e m w i t h the t r u e g e n e r a t i n g p r o c e s s .
Further , the e s t i m a t e s for the r a n d o m i n t e r c e p t s had to be c h a n g e d as w e l l to be a b l e to c o m p a r e t h e m to the t r u e v a l u e s .
It w i l l l a t e r be e x p l a i n e d how the t r a n s f o r m a t i o n is d o n e .
# In total , 2 0 8 . 0 0 0 M C M C s a m p l e s are drawn , w i t h a b u r n i n of 8 . 0 0 0 and a t h i n i n g of 400 , r e s u l t i n g in 500 s a m p l e s . A f t e r the M C M C a l g o r i t h m was d o n e w i t h BayesX , the d r a w s for e a c h p a r a m e t e r w e r e s a v e d to a f i l e . T h o s e f i l e s w i l l now be r e a d i n t o R and f u r t h e r a n a l y z e d .
#
-# R e a d in M C M C d r a w s
# Mu
mu . i n t e r c e p t < - r e a d.t a b l e(p a s t e 0(g e t w d() ," / o u t p u t _ M A I N _ mu _ R E G R E S S I O N _ yn _ L i n e a r E f f e c t s _ s a m p l e . raw ") , h e a der =T R U E) mu . s p l i n e < - r e a d.t a b l e(p a s t e 0(g e t w d() ," / o u t p u t _ M A I N _ mu _
R E G R E S S I O N _ yn _ n o n l i n e a r _ p s p l i n e _ e f f e c t _ of _ x _ s a m p l e . raw ") , h e a d er = T R U E)
mu . RI < - r e a d.t a b l e(p a s t e 0(g e t w d() ," / o u t p u t _ M A I N _ mu _ R E G R E S S I O N _ yn _ r a n d o m _ e f f e c t _ of _ g r o u p s _ s a m p l e . raw ") , h e a der = T R U E)
mu . var . RI < - r e a d.t a b l e(p a s t e 0(g e t w d() ," / o u t p u t _ M A I N _ mu _
R E G R E S S I O N _ yn _ v a r i a n c e _ of _ r a n d o m _ e f f e c t _ of _ g r o u p s _ s a m p l e . raw ") , h e a der = T R U E)
# S i g m a
s i g m a . i n t e r c e p t < - r e a d.t a b l e(p a s t e 0(g e t w d() ," / o u t p u t _ M A I N _ s i g m a _ R E G R E S S I O N _ yn _ L i n e a r E f f e c t s _ s a m p l e . raw ") , h e a der =T R U E)
s i g m a . s p l i n e < - r e a d.t a b l e(p a s t e 0(g e t w d() ," / o u t p u t _ M A I N _ s i g m a _ R E G R E S S I O N _ yn _ n o n l i n e a r _ p s p l i n e _ e f f e c t _ of _ x _ s a m p l e . raw ") , h e a d er = T R U E)
s i g m a . RI < - r e a d.t a b l e(p a s t e 0(g e t w d() ," / o u t p u t _ M A I N _ s i g m a _
R E G R E S S I O N _ yn _ r a n d o m _ e f f e c t _ of _ g r o u p s _ s a m p l e . raw ") , h e a der = T R U E)
s i g m a . var . RI < - r e a d.t a b l e(p a s t e 0(g e t w d() ," / o u t p u t _ M A I N _ s i g m a _ R E G R E S S I O N _ yn _ v a r i a n c e _ of _ r a n d o m _ e f f e c t _ of _ g r o u p s _ s a m p l e . raw ") , h e a der = T R U E)
# Nu
nu . i n t e r c e p t < - r e a d.t a b l e(p a s t e 0(g e t w d() ," / o u t p u t _ M A I N _ pi _ R E G R E S S I O N _ yn _ L i n e a r E f f e c t s _ s a m p l e . raw ") , h e a der =T R U E) nu . s p l i n e < - r e a d.t a b l e(p a s t e 0(g e t w d() ," / o u t p u t _ M A I N _ pi _
R E G R E S S I O N _ yn _ n o n l i n e a r _ p s p l i n e _ e f f e c t _ of _ x _ s a m p l e . raw ") , h e a d er = T R U E)
nu . RI < - r e a d.t a b l e(p a s t e 0(g e t w d() ," / o u t p u t _ M A I N _ pi _ R E G R E S S I O N _ yn _ r a n d o m _ e f f e c t _ of _ g r o u p s _ s a m p l e . raw ") , h e a der = T R U E)
nu . var . RI < - r e a d.t a b l e(p a s t e 0(g e t w d() ," / o u t p u t _ M A I N _ pi _
R E G R E S S I O N _ yn _ v a r i a n c e _ of _ r a n d o m _ e f f e c t _ of _ g r o u p s _ s a m p l e . raw ") , h e a der = T R U E)
#
-# 3.1 S a m p l i n g p a t h s and a u t o c o r r e l a t i o n s
# Firstly , the s a m p l i n g p a t h s and the a u t o c o r r e l a t i o n s of the d r a w s for any p a r a m e t e r w i l l be s h o w n .
# -f i l e s . a u t o < - l i s t . f i l e s( p a t t e r n = " \\_ a u t o c o r . raw $ ")
f i l e s . sam < - l i s t . f i l e s( p a t t e r n = " \\_ s a m p l e . raw $ ")
# L i n e a r E f f e c t s
# # S a m p l e P a t h s
f i l e s . l i n e a r e f f e c t < - f i l e s . sam [!is.na(s t r _ m a t c h( f i l e s . sam , "
L i n e a r E f f e c t s ") ) ]
lin .f r a m e < - m a t r i x(n c o l = 3 ,n r o w = 0) for( i in 1:l e n g t h( f i l e s . l i n e a r e f f e c t ) ) {
s a m p l e < - r e a d.t a b l e( f i l e s . l i n e a r e f f e c t [ i ] , h e a der = T R U E) p a r a m e t e r < - s t r _ m a t c h( f i l e s . l i n e a r e f f e c t [ i ] , " M A I N _ ( . * ? ) _ ") [1 ,2]
if( p a r a m e t e r == " pi ") p a r a m e t e r < - " nu "
t y p e < - s t r _ m a t c h( f i l e s . l i n e a r e f f e c t [ i ] , " yn _ ( . * ? ) _ ") [1 ,2]
if( p a r a m e t e r == " s i g m a ") {
f i r s t < - c b i n d( - s a m p l e [ ,2]/2 , type , p a r a m e t e r ) }e l s e{
if( p a r a m e t e r == " nu ") {
f i r s t < - c b i n d( - s a m p l e [ ,2] , type , p a r a m e t e r ) }e l s e{
f i r s t < - c b i n d( s a m p l e [ ,2] , type , p a r a m e t e r ) }
}
lin .f r a m e < - r b i n d( lin .frame, f i r s t ) }
lin .f r a m e < - d a t a.f r a m e( s a m p l e = as.n u m e r i c( lin .f r a m e[ ,1]) , T y p e = as.f a c t o r( lin .f r a m e[ ,2]) ,
P a r a m e t e r = f a c t o r( lin .f r a m e[ ,3] , l a b e l s = c(" mu ", " s i g m a ", " nu ") , l e v e l s = c(" mu ", " s i g m a ", " nu ") ) ,
x = 1:(n r o w( lin .f r a m e) /3) )
x y p l o t( s a m p l e ~ x | P a r a m e t e r , d a t a = lin .frame, l a y o u t = c(1 ,3) , par . s t r i p .t e x t = l i s t( cex = 1 . 2 ) , t y p e = " l ", m a i n = " S a m p l e
p a t h s - I n t e r c e p t s ", s c a l e s =l i s t( a l t e r n a t i n g = T , y =l i s t( r e l a t i o n =" f r e e ") , cex = 1 . 1 ) ,i n d e x. c o n d =l i s t(c(3 ,2 ,1) ) ,
p a r . s e t t i n g s=l i s t( par . m a i n .t e x t=l i s t( cex = 2 . 5 ) ) , x l a b = l i s t( l a b e l = " S a m p l e ID ", cex = 2) , y l a b = " ")
# # A u t o c o r r e l a t i o n
f i l e s . l i n e a r e f f e c t < - f i l e s . a u t o[!is.na(s t r _ m a t c h(f i l e s . a u t o, "
L i n e a r E f f e c t s ") ) ]
lin .f r a m e < - m a t r i x(n c o l = 3 ,n r o w = 0) for( i in 1:l e n g t h( f i l e s . l i n e a r e f f e c t ) ) {
a u t o c o r < - r e a d.t a b l e( f i l e s . l i n e a r e f f e c t [ i ] , h e a der = T R U E) p a r a m e t e r < - s t r _ m a t c h( f i l e s . l i n e a r e f f e c t [ i ] , " M A I N _ ( . * ? ) _ ") [1 ,2]
if( p a r a m e t e r == " pi ") p a r a m e t e r < - " nu "
t y p e < - s t r _ m a t c h( f i l e s . l i n e a r e f f e c t [ i ] , " yn _ ( . * ? ) _ ") [1 ,2]
f i r s t < - c b i n d(c(1 , a u t o c o r [ ,2]) , type , p a r a m e t e r ) lin .f r a m e < - r b i n d( lin .frame, f i r s t )
}
lin .f r a m e < - d a t a.f r a m e( A u t o c o r = as.n u m e r i c( lin .f r a m e[ ,1]) , T y p e = as.f a c t o r( lin .f r a m e[ ,2]) ,
P a r a m e t e r = f a c t o r( lin .f r a m e[ ,3] , l a b e l s = c(" mu ", " s i g m a ", " nu ") , l e v e l s = c(" mu ", " s i g m a ", " nu ") ) ,
x = 1 : 5 1 )
x y p l o t( A u t o c o r ~ x | P a r a m e t e r , d a t a = lin .frame, t y p e = " h ", l a y o u t
= c(3 ,1) , par . s t r i p .t e x t = l i s t( cex = 2) , m a i n = "
A u t o c o r r e l a t i o n - I n t e r c e p t s ", x l a b = l i s t( l a b e l = " Lag ", cex = 2) , y l a b = l i s t( l a b e l = " A u t o c o r r e l a t i o n ", cex = 2) , s c a l e s = l i s t( cex = 1 . 1 ) , p a r . s e t t i n g s=l i s t( par . m a i n .t e x t=l i s t( cex
= 2 . 5 ) ) )
# B a s i s C o e f f i c i e n t s
# # S a m p l e P a t h s
f i l e s . n o n l i n e a r < - f i l e s . sam [!is.na(s t r _ m a t c h( f i l e s . sam , " yn _ n o n l i n e a r ") ) ]
n o n l i n .f r a m e < - m a t r i x(n c o l = 4 ,n r o w = 0) for( i in 1:l e n g t h( f i l e s . n o n l i n e a r ) ) {
s a m p l e < - r e a d.t a b l e( f i l e s . n o n l i n e a r [ i ] , h e a der = T R U E)
p a r a m e t e r < - s t r _ m a t c h( f i l e s . n o n l i n e a r [ i ] , " M A I N _ ( . * ? ) _ ") [1 ,2]
if( p a r a m e t e r == " pi ") p a r a m e t e r < - " nu "
t y p e < - s t r _ m a t c h( f i l e s . n o n l i n e a r [ i ] , " yn _ ( . * ? ) _ ") [1 ,2]
for( j in 1 : 2 2 ) {
if( p a r a m e t e r == " s i g m a ") {
f i r s t < - c b i n d( - s a m p l e [ ,1+ j ]/2 , type , p a r a m e t e r , r a n d I D = j ) }e l s e{
if( p a r a m e t e r == " nu ") {
f i r s t < - c b i n d( - s a m p l e [ ,1+ j ] , type , p a r a m e t e r , r a n d I D = j ) }e l s e{
f i r s t < - c b i n d( s a m p l e [ ,1+ j ] , type , p a r a m e t e r , r a n d I D = j ) }
}
n o n l i n .f r a m e < - r b i n d( n o n l i n .frame, f i r s t ) }
}
n o n l i n .f r a m e < - d a t a.f r a m e( s a m p l e = as.n u m e r i c( n o n l i n .f r a m e[ ,1]) , T y p e = as.f a c t o r( n o n l i n .f r a m e[ ,2]) ,
P a r a m e t e r = as.f a c t o r( n o n l i n .f r a m e[ ,3]) ,
r a n d I D = f a c t o r( n o n l i n .f r a m e[ ,4] , l e v e l s = as.c h a r a c t e r( 1 : 2 2 ) , l a b e l s = as.c h a r a c t e r( 1 : 2 2 ) ) ,
x = 1 : 5 0 0 )
x y p l o t( s a m p l e ~ x | r a n d I D + P a r a m e t e r , d a t a = n o n l i n .f r a m e[ n o n l i n . f r a m e$ P a r a m e t e r == " mu ",] , t y p e = " l ", l a y o u t = c(4 ,6) , m a i n =
" S a m p l e p a t h s - B a s i s C o e f f i c i e n t s ", as.t a b l e=TRUE, x l a b = l i s t( l a b e l = " S a m p l e ID ", cex = 2) , y l a b = " ", p a r . s e t t i n g s=
l i s t( par . m a i n .t e x t=l i s t( cex = 2 . 5 ) ) , par . s t r i p .t e x t = l i s t( cex = 1 . 2 ) , s c a l e s = l i s t( cex = 1 . 1 ) )
x y p l o t( s a m p l e ~ x | r a n d I D + P a r a m e t e r , d a t a = n o n l i n .f r a m e[ n o n l i n . f r a m e$ P a r a m e t e r == " s i g m a ",] , t y p e = " l ", l a y o u t = c(4 ,6) , m a i n = " S a m p l e p a t h s - B a s i s C o e f f i c i e n t s ", as.t a b l e=TRUE, x l a b = l i s t( l a b e l = " S a m p l e ID ", cex = 2) , y l a b = " ",
p a r . s e t t i n g s=l i s t( par . m a i n .t e x t=l i s t( cex = 2 . 5 ) ) , par . s t r i p .t e x t
= l i s t( cex = 1 . 2 ) , s c a l e s = l i s t( cex = 1 . 1 ) )
x y p l o t( s a m p l e ~ x | r a n d I D + P a r a m e t e r , d a t a = n o n l i n .f r a m e[ n o n l i n . f r a m e$ P a r a m e t e r == " nu ",] , t y p e = " l ", l a y o u t = c(4 ,6) , m a i n =
" S a m p l e p a t h s - B a s i s C o e f f i c i e n t s ", as.t a b l e=TRUE, x l a b = l i s t( l a b e l = " S a m p l e ID ", cex = 2) , y l a b = " ", p a r . s e t t i n g s= l i s t( par . m a i n .t e x t=l i s t( cex = 2 . 5 ) ) , par . s t r i p .t e x t = l i s t( cex =
1 . 2 ) , s c a l e s = l i s t( cex = 1 . 1 ) )
# # A u t o c o r r e l a t i o n
f i l e s . r a n d o m < - f i l e s . a u t o[!is.na(s t r _ m a t c h(f i l e s . a u t o, " yn _ n o n l i n e a r ") ) ]
n o n l i n .f r a m e < - m a t r i x(n c o l = 4 ,n r o w = 0) for( i in 1:l e n g t h( f i l e s . r a n d o m ) ) {
a u t o c o r < - r e a d.t a b l e( f i l e s . r a n d o m [ i ] , h e a der = T R U E)
p a r a m e t e r < - s t r _ m a t c h( f i l e s . r a n d o m [ i ] , " M A I N _ ( . * ? ) _ ") [1 ,2]
if( p a r a m e t e r == " pi ") p a r a m e t e r < - " nu "
t y p e < - s t r _ m a t c h( f i l e s . r a n d o m [ i ] , " yn _ ( . * ? ) _ ") [1 ,2]
for( j in 1 : 2 2 ) {
f i r s t < - c b i n d(c(1 , a u t o c o r [ ,1+ j ]) , type , p a r a m e t e r , r a n d I D
= j )
n o n l i n .f r a m e < - r b i n d( n o n l i n .frame, f i r s t ) }
}
n o n l i n .f r a m e < - d a t a.f r a m e( A u t o c o r = as.n u m e r i c( n o n l i n .f r a m e[ ,1]) ,
T y p e = as.f a c t o r( n o n l i n .f r a m e[ ,2]) ,
P a r a m e t e r = as.f a c t o r( n o n l i n .f r a m e[ ,3]) ,
r a n d I D = f a c t o r( n o n l i n .f r a m e[ ,4] , l e v e l s = as.c h a r a c t e r( 1 : 2 2 ) , l a b e l s = as.c h a r a c t e r( 1 : 2 2 ) ) ,
x = 1 : 5 1 )
x y p l o t( A u t o c o r ~ x | r a n d I D + P a r a m e t e r , d a t a = n o n l i n .f r a m e[ n o n l i n . f r a m e$ P a r a m e t e r == " mu ",] , t y p e = " h ", l a y o u t = c(4 ,6) , m a i n =
" A u t o c o r r e l a t i o n - B a s i s C o e f f i c i e n t s ", as.t a b l e=TRUE, x l a b = l i s t( l a b e l = " Lag ", cex = 2) , y l a b = l i s t( l a b e l = "
A u t o c o r r e l a t i o n ", cex = 2) , p a r . s e t t i n g s=l i s t( par . m a i n .t e x t= l i s t( cex = 2 . 5 ) ) , par . s t r i p .t e x t = l i s t( cex = 1 . 2 ) , s c a l e s = l i s t( cex = 1 . 1 ) )
x y p l o t( A u t o c o r ~ x | r a n d I D + P a r a m e t e r , d a t a = n o n l i n .f r a m e[ n o n l i n . f r a m e$ P a r a m e t e r == " s i g m a ",] , t y p e = " h ", l a y o u t = c(4 ,6) , m a i n = " A u t o c o r r e l a t i o n - B a s i s C o e f f i c i e n t s ", as.t a b l e=TRUE, x l a b = l i s t( l a b e l = " Lag ", cex = 2) , y l a b = l i s t( l a b e l = "
A u t o c o r r e l a t i o n ", cex = 2) , p a r . s e t t i n g s=l i s t( par . m a i n .t e x t= l i s t( cex = 2 . 5 ) ) , par . s t r i p .t e x t = l i s t( cex = 1 . 2 ) , s c a l e s = l i s t( cex = 1 . 1 ) )
x y p l o t( A u t o c o r ~ x | r a n d I D + P a r a m e t e r , d a t a = n o n l i n .f r a m e[ n o n l i n . f r a m e$ P a r a m e t e r == " nu ",] , t y p e = " h ", l a y o u t = c(4 ,6) , m a i n =
" A u t o c o r r e l a t i o n - B a s i s C o e f f i c i e n t s ", as.t a b l e=TRUE, x l a b = l i s t( l a b e l = " Lag ", cex = 2) , y l a b = l i s t( l a b e l = "
A u t o c o r r e l a t i o n ", cex = 2) , p a r . s e t t i n g s=l i s t( par . m a i n .t e x t= l i s t( cex = 2 . 5 ) ) , par . s t r i p .t e x t = l i s t( cex = 1 . 2 ) , s c a l e s = l i s t( cex = 1 . 1 ) )
# R a n d o m I n t e r c e p t s
# # S a m p l e P a t h s
f i l e s . r a n d o m < - f i l e s . sam [!is.na(s t r _ m a t c h( f i l e s . sam , " yn _ r a n d o m "
) ) ]
ran .f r a m e < - m a t r i x(n c o l = 4 ,n r o w = 0) for( i in 1:l e n g t h( f i l e s . r a n d o m ) ) {
s a m p l e < - r e a d.t a b l e( f i l e s . r a n d o m [ i ] , h e a der = T R U E)
p a r a m e t e r < - s t r _ m a t c h( f i l e s . r a n d o m [ i ] , " M A I N _ ( . * ? ) _ ") [1 ,2]
if( p a r a m e t e r == " pi ") p a r a m e t e r < - " nu "
t y p e < - s t r _ m a t c h( f i l e s . r a n d o m [ i ] , " yn _ ( . * ? ) _ ") [1 ,2]
for( j in 1 : 1 0 0 ) {
if( p a r a m e t e r == " s i g m a ") {
f i r s t < - c b i n d( - s a m p l e [ ,1+ j ]/2 , type , p a r a m e t e r , r a n d I D = j ) }e l s e{
if( p a r a m e t e r == " nu ") {
f i r s t < - c b i n d( - s a m p l e [ ,1+ j ] , type , p a r a m e t e r , r a n d I D = j ) }e l s e{
f i r s t < - c b i n d( s a m p l e [ ,1+ j ] , type , p a r a m e t e r , r a n d I D = j ) }
}
ran .f r a m e < - r b i n d( ran .frame, f i r s t ) }
}
ran .f r a m e < - d a t a.f r a m e( s a m p l e = as.n u m e r i c( ran .f r a m e[ ,1]) , T y p e = as.f a c t o r( ran .f r a m e[ ,2]) ,
P a r a m e t e r = as.f a c t o r( ran .f r a m e[ ,3]) ,
r a n d I D = f a c t o r( ran .f r a m e[ ,4] , l e v e l s = as.c h a r a c t e r( 1 : 1 0 0 ) , l a b e l s = as.c h a r a c t e r( 1 : 1 0 0 ) ) ,
x = 1 : 5 0 0 )
x y p l o t( s a m p l e ~ x | r a n d I D + P a r a m e t e r , d a t a = ran .frame, t y p e = " l ", l a y o u t = c(10 ,10) , m a i n = " S a m p l e p a t h s - R a n d o m I n t e r c e p t s ", as.t a b l e=TRUE, y l i m = c( -3 ,3) , x l a b = l i s t( l a b e l = " S a m p l e ID "
, cex = 2) , y l a b = " ", par . s t r i p .t e x t = l i s t( cex = 1 . 2 ) ,
p a r . s e t t i n g s=l i s t( par . m a i n .t e x t=l i s t( cex = 2 . 5 ) ) , s c a l e s = l i s t( cex = 1 . 1 ) )
# # A u t o c o r r e l a t i o n
f i l e s . r a n d o m < - f i l e s . a u t o[!is.na(s t r _ m a t c h(f i l e s . a u t o, " yn _ r a n d o m ") ) ]
ran .f r a m e < - m a t r i x(n c o l = 4 ,n r o w = 0) for( i in 1:l e n g t h( f i l e s . r a n d o m ) ) {
a u t o c o r < - r e a d.t a b l e( f i l e s . r a n d o m [ i ] , h e a der = T R U E)
p a r a m e t e r < - s t r _ m a t c h( f i l e s . r a n d o m [ i ] , " M A I N _ ( . * ? ) _ ") [1 ,2]
if( p a r a m e t e r == " pi ") p a r a m e t e r < - " nu "
t y p e < - s t r _ m a t c h( f i l e s . r a n d o m [ i ] , " yn _ ( . * ? ) _ ") [1 ,2]
for( j in 1 : 1 0 0 ) {
f i r s t < - c b i n d(c(1 , a u t o c o r [ ,1+ j ]) , type , p a r a m e t e r , r a n d I D
= j )
ran .f r a m e < - r b i n d( ran .frame, f i r s t ) }
}
ran .f r a m e < - d a t a.f r a m e( A u t o c o r = as.n u m e r i c( ran .f r a m e[ ,1]) , T y p e = as.f a c t o r( ran .f r a m e[ ,2]) ,
P a r a m e t e r = as.f a c t o r( ran .f r a m e[ ,3]) ,
r a n d I D = f a c t o r( ran .f r a m e[ ,4] , l e v e l s = as.c h a r a c t e r( 1 : 1 0 0 ) , l a b e l s = as.c h a r a c t e r( 1 : 1 0 0 ) ) ,
x = 1 : 5 1 )
x y p l o t( A u t o c o r ~ x | r a n d I D + P a r a m e t e r , d a t a = ran .frame, t y p e = " h "
, l a y o u t = c(10 ,10) , m a i n = " A u t o c o r r e l a t i o n - R a n d o m I n t e r c e p t s ", as.t a b l e=TRUE, y l a b = l i s t( l a b e l = "
A u t o c o r r e l a t i o n ", cex = 2) , p a r . s e t t i n g s=l i s t( par . m a i n .t e x t= l i s t( cex = 2 . 5 ) ) , x l a b = l i s t( l a b e l = " Lag ", cex = 2) , par . s t r i p . t e x t = l i s t( cex = 1 . 2 ) , s c a l e s = l i s t( cex = 1 . 1 ) )
# B a s i s C o e f f i c i e n t s V a r i a n c e
# # S a m p l e P a t h s
f i l e s . v a r n o n l i n e a r < - f i l e s . sam [!is.na(s t r _ m a t c h( f i l e s . sam , " yn _ v a r i a n c e _ of _ n o n l i n e a r _ p s p l i n e ") ) ]
v a r n o n l i n .f r a m e < - m a t r i x(n c o l = 3 ,n r o w = 0) for( i in 1:l e n g t h( f i l e s . v a r n o n l i n e a r ) ) {
s a m p l e < - r e a d.t a b l e( f i l e s . v a r n o n l i n e a r [ i ] , h e a der = T R U E)
p a r a m e t e r < - s t r _ m a t c h( f i l e s . v a r n o n l i n e a r [ i ] , " M A I N _ ( . * ? ) _ ") [1 ,2]
if( p a r a m e t e r == " pi ") p a r a m e t e r < - " nu "
t y p e < - s t r _ m a t c h( f i l e s . v a r n o n l i n e a r [ i ] , " yn _ ( . * ? ) _ ") [1 ,2]
f i r s t < - c b i n d( s a m p l e [ ,2] , type , p a r a m e t e r )
v a r n o n l i n .f r a m e < - r b i n d( v a r n o n l i n .frame, f i r s t ) }
v a r n o n l i n .f r a m e < - d a t a.f r a m e( s a m p l e = as.n u m e r i c( v a r n o n l i n .f r a m e [ ,1]) ,
T y p e = as.f a c t o r( v a r n o n l i n .f r a m e[ ,2]) ,
P a r a m e t e r = f a c t o r( v a r n o n l i n .f r a m e[ ,3] , l a b e l s = c(" mu ", " s i g m a ",
" nu ") , l e v e l s = c(" mu ", " s i g m a ", " nu ") ) , x = 1:(n r o w( v a r n o n l i n .f r a m e) /3) )
x y p l o t( s a m p l e ~ x | P a r a m e t e r , d a t a = v a r n o n l i n .frame, l a y o u t = c (1 ,3) , t y p e = " l ", m a i n = " S a m p l e p a t h s - V a r i a n c e B a s i s C o e f f i c i e n t s ", s c a l e s =l i s t( a l t e r n a t i n g = T , cex = 1.1 , y =l i s t( r e l a t i o n =" f r e e ") ) ,i n d e x. c o n d =l i s t(c(3 ,2 ,1) ) , x l a b = l i s t( l a b e l
= " S a m p l e ID ", cex = 2) , y l a b = " ", p a r . s e t t i n g s=l i s t( par . m a i n .t e x t=l i s t( cex = 2 . 5 ) ) , par . s t r i p .t e x t = l i s t( cex = 1 . 2 ) )
# # A u t o c o r r e l a t i o n
f i l e s . n o n l i n e a r e f f e c t < - f i l e s . a u t o[!is.na(s t r _ m a t c h(f i l e s . a u t o,
" yn _ v a r i a n c e _ of _ n o n l i n e a r _ p s p l i n e ") ) ] v a r n o n l i n .f r a m e < - m a t r i x(n c o l = 3 ,n r o w = 0) for( i in 1:l e n g t h( f i l e s . n o n l i n e a r e f f e c t ) ) {
a u t o c o r < - r e a d.t a b l e( f i l e s . n o n l i n e a r e f f e c t [ i ] , h e a der = T R U E) p a r a m e t e r < - s t r _ m a t c h( f i l e s . n o n l i n e a r e f f e c t [ i ] , " M A I N _ ( . * ? ) _ "
) [1 ,2]
if( p a r a m e t e r == " pi ") p a r a m e t e r < - " nu "
t y p e < - s t r _ m a t c h( f i l e s . n o n l i n e a r e f f e c t [ i ] , " yn _ ( . * ? ) _ ") [1 ,2]
f i r s t < - c b i n d(c(1 , a u t o c o r [ ,2]) , type , p a r a m e t e r ) v a r n o n l i n .f r a m e < - r b i n d( v a r n o n l i n .frame, f i r s t ) }
v a r n o n l i n .f r a m e < - d a t a.f r a m e( A u t o c o r = as.n u m e r i c( v a r n o n l i n . f r a m e[ ,1]) ,
T y p e = as.f a c t o r( v a r n o n l i n .f r a m e[ ,2]) ,
P a r a m e t e r = f a c t o r( v a r n o n l i n .f r a m e[ ,3] , l a b e l s = c(" mu ", " s i g m a ",
" nu ") , l e v e l s = c(" mu ", " s i g m a ", " nu ") ) , x = 1 : 5 1 )
x y p l o t( A u t o c o r ~ x | P a r a m e t e r , d a t a = v a r n o n l i n .frame, t y p e = " h ", l a y o u t = c(3 ,1) , m a i n = " A u t o c o r r e l a t i o n - V a r i a n c e B a s i s C o e f f i c i e n t s ", x l a b = l i s t( l a b e l = " Lag ", cex = 2) , y l a b = l i s t( l a b e l = " A u t o c o r r e l a t i o n ", cex = 2) , p a r . s e t t i n g s=l i s t( par . m a i n .t e x t=l i s t( cex = 2 . 5 ) ) , par . s t r i p .t e x t = l i s t( cex = 1 . 2 ) , s c a l e s = l i s t( cex = 1 . 1 ) )
# V a r i a n c e R a n d o m I n t e r c e p t s
# # S a m p l e P a t h s
f i l e s . v a r r a n d o m < - f i l e s . sam [!is.na(s t r _ m a t c h( f i l e s . sam , " yn _ v a r i a n c e _ of _ r a n d o m _ e f f e c t ") ) ]
v a r r a n d o m .f r a m e < - m a t r i x(n c o l = 3 ,n r o w = 0) for( i in 1:l e n g t h( f i l e s . v a r r a n d o m ) ) {
s a m p l e < - r e a d.t a b l e( f i l e s . v a r r a n d o m [ i ] , h e a der = T R U E)
p a r a m e t e r < - s t r _ m a t c h( f i l e s . v a r r a n d o m [ i ] , " M A I N _ ( . * ? ) _ ") [1 ,2]
if( p a r a m e t e r == " pi ") p a r a m e t e r < - " nu "
t y p e < - s t r _ m a t c h( f i l e s . v a r r a n d o m [ i ] , " yn _ ( . * ? ) _ ") [1 ,2]
if( p a r a m e t e r == " s i g m a ") {
f i r s t < - c b i n d( s a m p l e [ ,3] , type , p a r a m e t e r ) p r i n t(m e a n(as.n u m e r i c( f i r s t [ ,1]) ) )
}e l s e{
f i r s t < - c b i n d( s a m p l e [ ,2] , type , p a r a m e t e r ) }
v a r r a n d o m .f r a m e < - r b i n d( v a r r a n d o m .frame, f i r s t ) }
v a r r a n d o m .f r a m e < - d a t a.f r a m e( s a m p l e = as.n u m e r i c( v a r r a n d o m .f r a m e [ ,1]) ,
T y p e = as.f a c t o r( v a r r a n d o m .f r a m e[ ,2]) ,
P a r a m e t e r = f a c t o r( v a r r a n d o m .f r a m e[ ,3] , l a b e l s = c(" mu ", " s i g m a ",
" nu ") , l e v e l s = c(" mu ", " s i g m a ", " nu ") ) , x = 1:(n r o w( v a r r a n d o m .f r a m e) /3) )
x y p l o t( s a m p l e ~ x | P a r a m e t e r , d a t a = v a r r a n d o m .frame, l a y o u t = c (1 ,3) , t y p e = " l ", m a i n = " S a m p l e p a t h s - V a r i a n c e R a n d o m I n t e r c e p t s ",
s c a l e s =l i s t( a l t e r n a t i n g = T , y =l i s t( r e l a t i o n =" f r e e ") , cex = 1 . 1 ) , i n d e x. c o n d =l i s t(c(3 ,2 ,1) ) , x l a b = l i s t( l a b e l = " S a m p l e ID ", cex = 2) , y l a b = " ", p a r . s e t t i n g s=l i s t( par . m a i n .t e x t=l i s t( cex
= 2 . 5 ) ) , par . s t r i p .t e x t = l i s t( cex = 1 . 2 ) )
# # A u t o c o r r e l a t i o n
f i l e s . v a r r a n d o m < - f i l e s . a u t o[!is.na(s t r _ m a t c h(f i l e s . a u t o, " yn _ v a r i a n c e _ of _ r a n d o m _ e f f e c t ") ) ]
v a r r a n d o m .f r a m e < - m a t r i x(n c o l = 3 ,n r o w = 0) for( i in 1:l e n g t h( f i l e s . v a r r a n d o m ) ) {
a u t o c o r < - r e a d.t a b l e( f i l e s . v a r r a n d o m [ i ] , h e a der = T R U E)
p a r a m e t e r < - s t r _ m a t c h( f i l e s . v a r r a n d o m [ i ] , " M A I N _ ( . * ? ) _ ") [1 ,2]
if( p a r a m e t e r == " pi ") p a r a m e t e r < - " nu "
t y p e < - s t r _ m a t c h( f i l e s . v a r r a n d o m [ i ] , " yn _ ( . * ? ) _ ") [1 ,2]
f i r s t < - c b i n d(c(1 , a u t o c o r [ ,2]) , type , p a r a m e t e r ) v a r r a n d o m .f r a m e < - r b i n d( v a r r a n d o m .frame, f i r s t ) }
v a r r a n d o m .f r a m e < - d a t a.f r a m e( A u t o c o r = as.n u m e r i c( v a r r a n d o m . f r a m e[ ,1]) ,
T y p e = as.f a c t o r( v a r r a n d o m .f r a m e[ ,2]) ,
P a r a m e t e r = f a c t o r( v a r r a n d o m .f r a m e[ ,3] , l a b e l s = c(" mu ", " s i g m a ",
" nu ") , l e v e l s = c(" mu ", " s i g m a ", " nu ") ) , x = 1 : 5 1 )
x y p l o t( A u t o c o r ~ x | P a r a m e t e r , d a t a = v a r r a n d o m .frame, t y p e = " h ", l a y o u t = c(3 ,1) , m a i n = " A u t o c o r r e l a t i o n - V a r i a n c e R a n d o m I n t e r c e p t s ", x l a b = l i s t( l a b e l = " Lag ", cex = 2) , y l a b = l i s t( l a b e l = " A u t o c o r r e l a t i o n ", cex = 2) , p a r . s e t t i n g s=l i s t( par . m a i n .t e x t=l i s t( cex = 2 . 5 ) ) , par . s t r i p .t e x t = l i s t( cex = 1 . 2 ) )
#
-# 3.2 M o d e l i n g the e f f e c t of x on the p a r a m e t e r s
#
-# G e n e r a l S p l i n e E s t i m a t e f u n c t i o n
e s t i m a t e . s p l i n e < - f u n c t i o n( x , s p l i n e . s a m p l e s ) {
# S e q u e n c e x
s e q u e n c e < - seq( min ( x ) , max ( x ) , l e n g t h. out = 2 0 0 0 )
# P o s t e r i o r M e a n of B a s i s C o e f f i c i e n t s p o s t .m e a n < - c o l M e a n s( s p l i n e . s a m p l e s )
# B - B a s i s m a t r i x
bs . mat < - bs . d e s i g n ( x = s e q u e n c e , xl = min ( x ) , xr = max ( x ) , s p l i n e . d e g r e e =3 , n b a s i s = 2 2 )
# P o s t e r i o r m e a n for s p l i n e
n o n l i n e a r .m e a n < - bs . mat %*% m a t r i x( p o s t .mean, n c o l = 1)
# P o i n t w i s e c r e d i b i l i t y i n t e r v a l s
s p l i n e .t < - bs . mat %*% t( s p l i n e . s a m p l e s )
n o n l i n e a r . q u a n t i l e s < - t(a p p l y( s p l i n e .t, 1 , f u n c t i o n( x ) q u a n t i l e( x , p r o b s = c( 0 . 0 2 5 , 0 . 9 7 5 ) ) ) )
r e t u r n(d a t a.f r a m e( x = s e q u e n c e , m e a n = n o n l i n e a r .mean, q 0 2 5 = n o n l i n e a r . q u a n t i l e s [ ,1] , q 9 7 5 = n o n l i n e a r . q u a n t i l e s [ ,2]) ) }
#
-# 3 . 2 . 1 M o d e l i n g r e s u l t s for mu
# The p a r a m e t r i z a t i o n d o e s not c h a n g e for the p a r a m e t e r mu . As a log l i n k f u n c t i o n is u s e d for the p a r a m e t e r , the d r a w s are e x p o n e n t i a t e d and t h e n the e s t i m a t e is c a l c u l a t e d .
# -mu . s a m p l e s < - exp( mu . s p l i n e [ ,2:23] + mu . i n t e r c e p t $ b _1)
mu . est < - e s t i m a t e . s p l i n e ( df $ x , mu . s a m p l e s )
p l o t( mu . est $ x , mu . est $mean, t y p e = " l ", y l i m = c(5 ,60) ,
m a i n = e x p r e s s i o n(" E f f e c t of x on y for " * mu ) , y l a b = " ", x l a b
= " x ")
l i n e s( mu . est $ x , mu . est $ q025 , lty =2) l i n e s( mu . est $ x , mu . est $ q975 , lty =2)
l i n e s( mu . est $ x , mu . l i n k ( mu . est $ x ,0) , col = " red ", lwd = 2)
#
-# 3 . 2 . 2 M o d e l i n g r e s u l t s for s i g m a
# As a l r e a d y m e n t i o n e d , B a y e s X u s e s a d i f f e r e n t p a r a m e t r i z a t i o n for the s e c o n d p a r a m e t e r . The c o n n e c t i o n b e t w e e n s i g m a t h a t was u s e d in the DGP and tau is tau = 1/ s i g m a ^2. Further , s i m i l a r l y to the DGP of sigma , a log f u n c t i o n is a p p l i e d for tau in B a y e s X . Hence , the f o l l o w i n g t r a n s f o r m a t i o n is u s e d :
# s i g m a ( x _{ i , j }; b o l d s y m b o l { b e t a }_2 , g a m m a _{2 , i }) ^2 = 1/ s q r t ( tau ( x _{ i , j }; b o l d s y m b o l { b e t a }_{ tau } , g a m m a _{ tau , i }) ) = 1/ s q r t ( exp ( x _{ i , j } ’ b o l d s y m b o l { b e t a }_{ tau } + g a m m a _{ tau , i }) )
# -s i g m a . -s a m p l e -s < - 1/s q r t(exp( s i g m a . i n t e r c e p t $ b _1 + s i g m a . s p l i n e
[ , 2 : 2 3 ] ) )
s i g m a . est < - e s t i m a t e . s p l i n e ( df $ x , s i g m a . s a m p l e s )
p l o t( s i g m a . est $ x , s i g m a . est $mean, t y p e = " l ", y l i m = c(0 ,0.5) , m a i n = e x p r e s s i o n(" E f f e c t of x on y for " * s i g m a ) , y l a b = " ",
x l a b = " x ")
l i n e s( s i g m a . est $ x , s i g m a . est $ q025 , lty =2) l i n e s( s i g m a . est $ x , s i g m a . est $ q975 , lty =2)
l i n e s( s i g m a . est $ x , s i g m a . l i n k ( s i g m a . est $ x ,0) , col = " red ")
#
-# 3 . 2 . 3 M o d e l i n g r e s u l t s for nu
# For nu , a d i f f e r e n t p a r a m e t r i z a t i o n is u s e d as w e l l . B a y e s X u s e s pi = 1 - nu and thus , the f o l l o w i n g t r a n s f o r m a t i o n is a p p l i e d to o b t a i n the e s t i m a t i o n on nu :
# nu ( x _{ i , j }; b o l d s y m b o l { b e t a }_3 , g a m m a _{3 , i }) = 1 - pi ( x _{ i , j };
b o l d s y m b o l { b e t a }_{ pi } , g a m m a _{ pi , i }) = 1 - f r a c { 1 } { 1 + exp ( - x _{ i , j } ’ b o l d s y m b o l { b e t a }_{ pi } + g a m m a _{ pi , i }) } = f r a c { 1 } { 1 + exp ( x _{ i , j } ’ b o l d s y m b o l { b e t a }_{ pi } + g a m m a _{ pi , i }) }
# -nu . s a m p l e s < - 1 / ( 1 +exp(( nu . i n t e r c e p t [ ,2] + nu . s p l i n e [ , 2 : 2 3 ] ) ) ) nu . est < - e s t i m a t e . s p l i n e ( df $ x , nu . s a m p l e s )
p l o t( nu . est $ x , nu . est $mean, t y p e = " l ", y l i m = c(0 ,0.6) ,
m a i n = e x p r e s s i o n(" E f f e c t of x on y for " * nu ) , y l a b = " ", x l a b
= " x ")
l i n e s( nu . est $ x , nu . est $ q025 , lty =2) l i n e s( nu . est $ x , nu . est $ q975 , lty =2)
l i n e s( nu . est $ x , nu . l i n k ( nu . est $ x ,0) , col = " red ")
#
-# 3.3 M o d e l i n g r e s u l t s for the r a n d o m i n t e r c e p t s
# The re - p a r a m e t r i z a t i o n a l s o n e e d s to be t a k e n i n t o a c c o u n t for the e s t i m a t i o n of the r a n d o m i n t e r c e p t s .
#
#
-# 3 . 3 . 1 R a n d o m I n t e r c e p t s of mu
# A g a i n n o t h i n g c h a n g e s for mu :
# -mu . RI . est < - c o l M e a n s( mu . RI [ ,2:( g +1) ])
mu . RI . est < - mu . RI . est -m e a n( mu . RI . est ) mat < - m a t r i x(c( RI [ ,1] , mu . RI . est ) , b y r o w = TRUE, n r o w = 2)
r o w n a m e s( mat ) < - c(" T r u e _ RI ", " E s t i m a t e d _ RI ") m e a n(a p p l y( mat , 2 , d i f f ) ^2)
#
-# 3 . 3 . 2 R a n d o m I n t e r c e p t s of s i g m a
# For sigma , the r a n d o m i n t e r c e p t s n e e d to be re - c a l c u l a t e by :
# " s i g m a ^2" = & "1/ tau "
# exp ( g a m m a _{2 , i }) ^2 = 1/ exp ( g a m m a _{ tau , i })
# exp ( g a m m a _{2 , i }) ^2 = exp ( - g a m m a _{ tau , i })
# g a m m a _{2 , i } = - f r a c { g a m m a _{ tau , i } } { 2 }
# -s i g m a . RI . e-st < - c o l M e a n s( - s i g m a . RI [ ,2:( g +1) ] / 2 )
s i g m a . RI . est < - ( s i g m a . RI . est -m e a n( s i g m a . RI . est ) ) mat < - m a t r i x(c( RI [ ,2] , s i g m a . RI . est ) ,
b y r o w = TRUE, n r o w = 2)
r o w n a m e s( mat ) < - c(" T r u e _ RI ", " E s t i m a t e d _ RI ") mat
m e a n(a p p l y( mat , 2 , d i f f ) ^2)
#
-# 3 . 3 . 3 R a n d o m I n t e r c e p t s of nu
# For nu , the f o l l o w i n g t r a n s f o r m t a i o n n e e d s to be d o n e :
# " nu = 1 - pi "
# f r a c { 1 } { 1 + exp ( - g a m m a _{3 , i }) }=1 - f r a c { 1 } { 1 + exp ( - g a m m a _{ pi , i }) }
# f r a c { 1 } { 1 + exp ( - g a m m a _{3 , i }) }= f r a c { 1 } { 1 + exp ( g a m m a _{ pi , i }) }
# g a m m a _{3 , i } = - g a m m a _{ pi , i }
# -nu . RI . est < - c o l M e a n s( - nu . RI [ ,2:( g +1) ])
nu . RI . est < - nu . RI . est -m e a n( nu . RI . est ) mat < - m a t r i x(c( RI [ ,3] , nu . RI . est ) , b y r o w = TRUE, n r o w = 2)
r o w n a m e s( mat ) < - c(" T r u e _ RI ", " E s t i m a t e d _ RI ") mat
m e a n(a p p l y( mat , 2 , d i f f ) ^2)
#
-# 3 . 3 . 4 M o d e l i n g r e s u l t s for the v a r i a n c e s of the r a n d o m i n t e r c e p t s
# -d a t a.f r a m e( mu = c(0.1 , c o l M e a n s( mu . var . RI [ ,2:3]) [ 1 ] ) ,
s i m g a = c(0.5 , c o l M e a n s( s i g m a . var . RI [ ,2:3]) [ 2 ] ) ,
nu = c(0.25 , c o l M e a n s( nu . var . RI [ ,2:3]) [ 1 ] ) , row.n a m e s = c(" T r u e ",
" E s t i m a t e s ") )
#
-# 3 . 3 . 5 V i s u a l i z a t i o n of the t r u e a g a i n s t the e s t i m a t e d r a n d o m i n t e r c e p t s
# -n u m b e r < - df % >%
g r o u p _ b y( g r o u p s ) % >%
s u m m a r i s e( n u m b e r = sum( y == 0) )
par ( m f r o w = c(1 ,4) , mar = c(7.1 , 9.1 , 8.1 , 6 . 1 ) , mgp =c(5 ,2 ,0) ) x l i m .1. min < - min (c( RI [ ,1] , mu . RI . est ) )
x l i m .1. max < - max (c( RI [ ,1] , mu . RI . est ) ) x l i m .2. min < - min (c( RI [ ,2] , s i g m a . RI . est ) ) x l i m .2. max < - max (c( RI [ ,2] , s i g m a . RI . est ) )
x l i m .3. min < - min (c( RI [ ,3] , nu . RI . est ) ) x l i m .3. max < - max (c( RI [ ,3] , nu . RI . est ) )
p l o t( mu . RI . est , RI [ ,1] , x l i m = c( x l i m .1. min , x l i m .1. max ) , y l i m = c( x l i m .1. min , x l i m .1. max ) , pch = 16 ,
x l a b = " E s t i m a t e d r a n d o m i n t e r c e p t s ", y l a b = " T r u e r a n d o m
i n t e r c e p t s ", cex . lab = 2.5 , x a x t = ’ n ’, y a x t = ’ n ’, m a i n = " ", cex = 2)
a x i s(1 ,c e x . a x i s= 2 . 5 ) a x i s(2 ,c e x . a x i s= 2 . 5 )
t i t l e(e x p r e s s i o n(a t o p(p a s t e(" E s t i m a t e d vs . t r u e r a n d o m i n t e r c e p t s for ", mu ) ) ) , cex . m a i n = 3)
a b l i n e( b = 1 , a = 0)
p l o t( s i g m a . RI . est , RI [ ,2] , x l i m = c( x l i m .2. min , x l i m .2. max ) , y l i m
= c( x l i m .2. min , x l i m .2. max ) , pch = 16 ,
x l a b = " E s t i m a t e d r a n d o m i n t e r c e p t s ", y l a b = " T r u e r a n d o m
i n t e r c e p t s ", cex . lab = 2.5 , x a x t = ’ n ’, y a x t = ’ n ’, m a i n = " ", cex = 2)
a x i s(1 ,c e x . a x i s= 2 . 5 ) a x i s(2 ,c e x . a x i s= 2 . 5 )
t i t l e(e x p r e s s i o n(a t o p(p a s t e(" E s t i m a t e d vs . t r u e r a n d o m i n t e r c e p t s for ", s i g m a ) ) ) , cex . m a i n = 3)
a b l i n e( b = 1 , a = 0)
max . col < - max ( n u m b e r $ n u m b e r )
p l o t( nu . RI . est , RI [ ,3] , x l i m = c( x l i m .3. min , x l i m .3. max ) , y l i m = c( x l i m .3. min , x l i m .3. max ) , pch = 16 , col = c(r a i n b o w( max . col ) ) [ n u m b e r $ n u m b e r +1] , x l a b = " E s t i m a t e d r a n d o m i n t e r c e p t s ", y l a b =
" T r u e r a n d o m i n t e r c e p t s ", cex . lab = 2.5 , x a x t = ’ n ’, y a x t = ’ n
’, m a i n = " ", cex = 2) a x i s(1 ,c e x . a x i s= 2 . 5 ) a x i s(2 ,c e x . a x i s= 2 . 5 )
t i t l e(e x p r e s s i o n(a t o p(p a s t e(" E s t i m a t e d vs . t r u e r a n d o m i n t e r c e p t s for ", nu ) ) ) , cex . m a i n = 3)
a b l i n e( b = 1 , a = 0) p l o t.new()
l e g e n d(" c e n t e r ", l e g e n d = p a s t e(0: max . col , " t i m e s of ", " z e r o ") , col = r a i n b o w( max . col ) , pch = 16 , cex = 3)
#
-# 3 . 3 . 6 C o r r e l a t i o n p l o t m a t r i x
# -ri < - c b i n d( ‘ RI mu ‘ = mu . RI . est , ‘ RI sigma ‘ = s i g m a . RI . est , ‘ RI
nu ‘ = nu . RI . est )
# C o r r e l a t i o n p l o t
l a y o u t(r b i n d(c(1 , rep(2 ,3) , rep(3 ,3) , rep(4 ,3) ) , c(5 , rep(6 ,3) , rep(7 ,3) , rep(8 ,3) ) ,
c(5 , rep(6 ,3) , rep(7 ,3) , rep(8 ,3) ) , c(5 , rep(6 ,3) , rep(7 ,3) , rep(8 ,3) ) , c(9 , rep(10 ,3) , rep(11 ,3) , rep(12 ,3) ) , c(9 , rep(10 ,3) , rep(11 ,3) , rep(12 ,3) ) , c(9 , rep(10 ,3) , rep(11 ,3) , rep(12 ,3) ) , c(13 , rep(14 ,3) , rep(15 ,3) , rep(16 ,3) ) , c(13 , rep(14 ,3) , rep(15 ,3) , rep(16 ,3) ) , c(13 , rep(14 ,3) , rep(15 ,3) , rep(16 ,3) ) , rep(17 ,10) ) )
par ( mar = c( 0 . 2 , 0 . 2 , 0 . 2 , 0 . 2 ) )
# 1. N o t h i n g
p l o t(c(0 , 1) , c(0 , 1) , ann = F , bty = ’ n ’, t y p e = ’ n ’, x a x t = ’ n ’ , y a x t = ’ n ’)
# 2. Mu
p l o t(c(0 , 1) , c(0 , 1) , ann = F , bty = ’ n ’, t y p e = ’ n ’, x a x t = ’ n ’ , y a x t = ’ n ’)
t e x t( x = 0.5 , y = 0.5 , e x p r e s s i o n( mu ) , cex = 5 , col = " b l a c k ")
# 3. S i g m a
p l o t(c(0 , 1) , c(0 , 1) , ann = F , bty = ’ n ’, t y p e = ’ n ’, x a x t = ’ n ’ , y a x t = ’ n ’)
t e x t( x = 0.5 , y = 0.5 , e x p r e s s i o n( s i g m a ) , cex = 5 , col = " b l a c k ")
# 4. Nu
p l o t(c(0 , 1) , c(0 , 1) , ann = F , bty = ’ n ’, t y p e = ’ n ’, x a x t = ’ n ’ , y a x t = ’ n ’)
t e x t( x = 0.5 , y = 0.5 , e x p r e s s i o n( nu ) , cex = 5 , col = " b l a c k ")
# 5. Mu
p l o t(c(0 , 1) , c(0 , 1) , ann = F , bty = ’ n ’, t y p e = ’ n ’, x a x t = ’ n ’ , y a x t = ’ n ’)
t e x t( x = 0.5 , y = 0.5 , e x p r e s s i o n( mu ) , cex = 5 , col = " b l a c k ")
# P a r a m e t e r
x < - seq( -5 ,5 ,0.01)
den .1 < - d e n s i t y( ri [ ,1]) den .2 < - d e n s i t y( ri [ ,2]) den .3 < - d e n s i t y( ri [ ,3])
x l i m .1 < - c( min ( ri [ ,1]) -0.2 , max ( ri [ ,1]) + 0 . 2 )
y l i m .1 < - c( min (c( den .1$ y ,d n o r m( x , sd = s q r t( 0 . 1 ) ) ) ) , max (c( den .1$ y ,d n o r m( x , sd = s q r t( 0 . 1 ) ) ) ) )
x l i m .2 < - c( min ( ri [ ,2]) -0.2 , max ( ri [ ,2]) + 0 . 2 )
y l i m .2 < - c( min (c( den .2$ y ,d n o r m( x , sd = s q r t( 0 . 5 ) ) ) ) , max (c( den .2$ y ,d n o r m( x , sd = s q r t( 0 . 5 ) ) ) ) )
x l i m .3 < - c( min ( ri [ ,3]) -0.2 , max ( ri [ ,3]) + 0 . 2 )
y l i m .3 < - c( min (c( den .3$ y , d n o r m( x , sd = s q r t( 0 . 2 5 ) ) ) ) , max (c( den .3$ y ,d n o r m( x , sd = s q r t( 0 . 2 5 ) ) ) ) )
# 6. RI mu
p l o t(d e n s i t y( ri [ ,1]) , x l i m = x l i m .1 , x a x t = ’ n ’, m a i n = " ", y a x t
= ’ n ’, lwd = 3 ,
y l i m = y l i m .1 , x l a b = " ", y l a b = " ") a x i s(3 , c e x . a x i s=2)
rug( ri [ ,1])
l i n e s( x , d n o r m( x , sd = s q r t( 0 . 1 ) ) , col = " b l u e ", lty = 2 , lwd = 3)
# par ( mar = c (0.2 ,0.2 ,4 ,0.2) )
# 7. Cor mu - s i g m a
p l o t( NA , x l i m = x l i m .1 , y l i m = y l i m .1 , y l a b = " ", x l a b = " ", x a x t = ’ n ’, y a x t = ’ n ’)
t e x t(m e a n( x l i m .1) , m e a n( y l i m .1) , l a b e l s = r o u n d(cor( ri ) [1 ,2] ,3) , cex = 5)
# 8. Cor mu - nu
p l o t( NA , x l i m = x l i m .2 , y l i m = y l i m .2 , y l a b = " ", x l a b = " ", x a x t = ’ n ’, y a x t = ’ n ’)
t e x t(m e a n( x l i m .2) , m e a n( y l i m .2) , l a b e l s = r o u n d(cor( ri ) [1 ,3] ,3) , cex = 5)
# 9. S i g m a
p l o t(c(0 , 1) , c(0 , 1) , ann = F , bty = ’ n ’, t y p e = ’ n ’, x a x t = ’ n ’ , y a x t = ’ n ’)
t e x t( x = 0.5 , y = 0.5 , e x p r e s s i o n( s i g m a ) , cex = 5 , col = " b l a c k ")
# 10. Dot mu s i g m a
c o e f.1 < - c o e f(lm( ri [ ,2]~ ri [ ,1]) )
p l o t( ri [ ,1] , ri [ ,2] , pch = 16 , x l a b = " ", y l a b = " ", x a x t = ’ n ’, cex = 1.5 ,c e x . a x i s=2 ,
x l i m = x l i m .1 , y l i m = x l i m .2)
a b l i n e( a = c o e f.1[1] , b = c o e f.1[2] , col = " red ", lwd = 2)
# 11. RI s i g m a
p l o t(d e n s i t y( ri [ ,2]) , x l i m = x l i m .2 , x a x t = ’ n ’, m a i n = " ", x a x t
= ’ n ’, y a x t = ’ n ’, lwd = 3 ,
y l i m = y l i m .2 , x l a b = " ", y l a b = " ") rug( ri [ ,2])
l i n e s( x , d n o r m( x , sd = s q r t( 0 . 5 ) ) , col = " b l u e ", lty = 2 , lwd = 3)
# 12. Cor sigma - nu
p l o t( NA , x l i m = x l i m .2 , y l i m = y l i m .2 , y l a b = " ", x l a b = " ", x a x t = ’ n ’, y a x t = ’ n ’)
t e x t(m e a n( x l i m .2) , m e a n( y l i m .2) , l a b e l s = r o u n d(cor( ri ) [2 ,3] ,3) , cex = 5)
# 13. Nu
p l o t(c(0 , 1) , c(0 , 1) , ann = F , bty = ’ n ’, t y p e = ’ n ’, x a x t = ’ n ’ , y a x t = ’ n ’)
t e x t( x = 0.5 , y = 0.5 , e x p r e s s i o n( nu ) , cex = 5 , col = " b l a c k ")
# 14. Dot mu nu
c o e f.2 < - c o e f(lm( ri [ ,3]~ ri [ ,1]) )
p l o t( ri [ ,1] , ri [ ,3] , pch = 16 , x l a b = " ", y l a b = " ", cex = 1.5 , c e x . a x i s=2 ,
x l i m = x l i m .1 , y l i m = x l i m .3)
a b l i n e( a = c o e f.2[1] , b = c o e f.2[2] , col = " red ", lwd = 2)
# 15. Dot s i g m a nu
c o e f.3 < - c o e f(lm( ri [ ,3]~ ri [ ,2]) )
p l o t( ri [ ,2] , ri [ ,3] , pch = 16 , x l a b = " ", y l a b = " ", y a x t = ’ n ’, cex = 1.5 ,c e x . a x i s=2 ,
x l i m = x l i m .2 , y l i m = x l i m .3)
a b l i n e( a = c o e f.3[1] , b = c o e f.3[2] , col = " red ", lwd = 2)
# 16. RI nu
p l o t(d e n s i t y( ri [ ,3]) , x l i m = x l i m .3 , x a x t = ’ n ’, m a i n = " ", y a x t
= ’ n ’, lwd = 3 ,c e x . a x i s=2 ,
y l i m = y l i m .3 , x l a b = " ", y l a b = " ") a x i s(1 ,c e x . a x i s=2)
rug( ri [ ,3])
l i n e s( x , d n o r m( x , sd = s q r t( 0 . 2 5 ) ) , col = " b l u e ", lty = 2 , lwd = 3)
par ( m f r o w = c(1 ,1) , mar = c(5.1 , 4.1 , 4.1 , 4 . 1 ) )
#
-# 3.4 F u r t h e r r e s u l t s
#
-# D e v i a n c e I n f o r m a t i o n C r i t e r i o n ( DIC )
# The DIC i t s e l f d o e s not t e l l us a n y t h i n g a b o u t the m o d e l but can be u s e d to c o m p a r e the m o d e l fit b e t w e e n m o d e l s ( the smaller , the b e t t e r ) .
dic < - r e a d.t a b l e(p a s t e 0(g e t w d() ," / o u t p u t _ M A I N _ mu _ R E G R E S S I O N _ y _ p r e d i c t _ DIC . res ") , h e a der = T R U E)
dic
# P r e d i c t i o n s of the f i r s t 5 o b s e r v a t i o n s :
p r e d < - r e a d.t a b l e(p a s t e 0(g e t w d() ," / o u t p u t _ M A I N _ mu _ R E G R E S S I O N _ y _ p r e d i c t . res ") , h e a der = T R U E)
h e a d( p r e d )
Bibliography
Belitz, C., Brezger, A., Klein, N., Kneib, T., Lang, S., and Umlauf, N. (2015a). BayesX -Software for Bayesian inference in structured additive regression models.
Belitz, C., Brezger, A., Klein, N., Kneib, T., Lang, S., and Umlauf, N. (2015b). BayesX - Software for Bayesian inference in structured additive regression models: Methodology Manual.
Belitz, C., Brezger, A., Klein, N., Kneib, T., Lang, S., and Umlauf, N. (2015c). BayesX -Software for Bayesian inference in structured additive regression models: Reference Man-ual.
Belitz, C., Brezger, A., Klein, N., Kneib, T., Lang, S., and Umlauf, N. (2015d). BayesX -Software for Bayesian inference in structured additive regression models: Tutorials.
Bishop, C. M. (2006). Pattern recognition and machine learning. Information science and statistics. Springer, New York.
Brezger, A. and Lang, S. (2006). Generalized structured additive regression based on bayesian p-splines. Comput. Statist. Data Anal, 50:967–991.
Casella, G. and Berger, R. L. (2002). Statistical inference. Duxbury advanced series. Brook-s/Cole, Cengage Learning, Belmont, Calif. OCLC: 935033840.
Dahl, D. B., Scott, D., Roosen, C., Magnusson, A., and Swinton, J. (2019). xtable: Export Tables to LaTeX or HTML. R package version 1.8-4.
Dowle, M. and Srinivasan, A. (2019). data.table: Extension of ‘data.frame‘. R package version 1.12.2.
Eilers, P. H. C. and Marx, B. D. (1996). Flexible smoothing with B -splines and penalties.
Statistical Science, 11(2):89–121.
Fahrmeir, L., Kneib, T., and Lang, S. (2004). Penalized structured additive regression for spacetime data: a bayesian perspective. Statistica Sinica, pages 731–761.
Fahrmeir, L., Kneib, T., and Lang, S. (2007). Regression: Modelle, Methoden und Anwen-dungen. Statistik und ihre AnwenAnwen-dungen. Springer, Berlin. OCLC: 180908931.
F¨ockersperger, T. and Wagner, H. (2020). Modeling mothers’ yearly earnings after returning from maternity leave with a bayesian distributional regression model. Proceedings of the 35th International Workshop on Statistical Modelling, pages 333–336.
Gamerman, D. (1997). Sampling from the posterior distribution in generalized linear mixed models. Statistics and Computing, 7:57–68.
Gelfand, A. E. and Smith, A. F. M. (1990). Sampling-based approaches to calculating marginal densities. Journal of the American Statistical Association, 85(410):398–409.
Gelman, A. (2005). Analysis of variance—why it is more important than ever. The Annals of Statistics, 33(1):1–53.
Genz, A., Bretz, F., Miwa, T., Mi, X., Leisch, F., Scheipl, F., and Hothorn, T. (2019).
mvtnorm: Multivariate Normal and t Distributions. R package version 1.0-10.
Grolemund, G. and Wickham, H. (2011). Dates and times made easy with lubridate. Journal of Statistical Software, 40(3):1–25.
Groll, A. (2016). PenCoxFrail: Regularization in Cox Frailty Models. R package version 1.0.1.
Hastie, T. J. and Tibshirani, R. J. (1990). Generalized additive models. London: Chapman
& Hall.
Hoff, P. D. (2009). A first course in Bayesian statistical methods. Springer texts in statistics.
Springer, London ; New York. OCLC: ocn310401109.
Jacobi, L., Wagner, H., and Fr¨uhwirth-Schnatter, S. (2016). Bayesian treatment effects models with variable selection for panel outcomes with an application to earnings effects of maternity leave. Journal of Econometrics, 193(1):234–250.
Jiang, J. (2007).Linear and Generalized Linear Mixed Models and Their Applications. Springer series in statistics. Springer, Dordrecht.
Kammann, E. E. and Wand, M. P. (2003). Geoadditive models. Journal of the Royal Statistical Society Series C, 52(1):1–18.
Klein, N. and Kneib, T. (2016). Scale-Dependent Priors for Variance Parameters in Struc-tured Additive Distributional Regression. Bayesian Analysis, 11(4):1071–1106.
Klein, N., Kneib, T., and Lang, S. (2013). Bayesian structured additive distributional re-gression. Working papers, Faculty of Economics and Statistics, University of Innsbruck.
Klein, N., Kneib, T., Lang, S., and Sohn, A. (2015a). Bayesian structured additive distri-butional regression with an application to regional income inequality in Germany. The Annals of Applied Statistics, 9(2):1024–1052.
Klein, N., Kneib, T., Lang, S., and Sohn, A. (2015b). Supplement to “Bayesian struc-tured additive distributional regression with an application to regional income inequality in Germany”.
Kneib, T., Klein, N., Lang, S., and Umlauf, N. (2019). Modular regression - a Lego sys-tem for building structured additive distributional regression models with tensor product interactions. TEST, 28(1):1–39.
Lalive, R. and Zweim¨uller, J. (2009). How does Parental Leave Affect Fertility and Return to Work? Evidence from Two Natural Experiments. Quarterly Journal of Economics, 124(3):1363–1402.
Lang, S. and Brezger, A. (2004). Bayesian P-Splines.Journal of Computational and Graphical Statistics, 13(1):183–212.
Lin, X. and Zhang, D. (1999). Inference in generalized additive mixed modelsby using smoothing splines. Journal of the Royal Statistical Society Series B, 61(2):381–400.
Metropolis, N., Rosenbluth, A. W., Rosenbluth, M. N., Teller, A. H., and Teller, E. (1953).
Equation of state calculations by fast computing machines. The Journal of Chemical Physics, 21(6):1087–1092.
Nelder, J. A. and Wedderburn, R. W. M. (1972). Generalized linear models. Journal of the Royal Statistical Society, Series A, General, 135:370–384.
R Core Team (2019). R: A Language and Environment for Statistical Computing. R Foun-dation for Statistical Computing, Vienna. OCLC: 1120300286.
Rigby, R. A. and Stasinopoulos, D. M. (2005). Generalized additive models for location, scale and shape (with discussion).Journal of the Royal Statistical Society: Series C (Applied Statistics), 54(3):507–554.
Robert, C. P. and Casella, G. (2010). Introducing Monte Carlo methods with R. Use R!
Springer, New York. OCLC: ocn462920377.
RStudio Team (2018). RStudio: Integrated Development Environment for R. RStudio, Inc., Boston, MA.
Sarkar, D. (2008). Lattice: Multivariate Data Visualization with R. Springer, New York. ISBN 978-0-387-75968-5.
Schmidt, P. (2017).Bayesian inference for structured additive regression models for large-scale problems with applications to medical imaging. PhD thesis, Ludwig-Maximilians-Universit¨at M¨unchen. Medium: application/pdf.
Signorell, A. et al. (2019). DescTools: Tools for Descriptive Statistics. R package version 0.99.28.
Singer, J. D. (1996). Applied longitudinal data analysis : modeling change and event occur-rence. Oxford U.P.
Umlauf, N., Adler, D., Kneib, T., Lang, S., and Zeileis, A. (2015). Structured Additive Regression Models: An R Interface to BayesX. Journal of Statistical Software, 63(21).
Umlauf, N., Klein, N., Simon, T., and Zeileis, A. (2019). bamlss: A Lego toolbox for flexible Bayesian regression (and beyond). arXiv 1909.12345, arXiv.org E-Print Archive.
Umlauf, N., Klein, N., and Zeileis, A. (2018). BAMLSS: Bayesian Additive Models for Loca-tion, Scale, and Shape (and Beyond). Journal of Computational and Graphical Statistics, 27(3):612–627.
Umlauf, N. and Kneib, T. (2018). A primer on Bayesian distributional regression. Statistical Modelling, 18(3-4):219–247.
Wickham, H. (2011). The split-apply-combine strategy for data analysis. Journal of Statis-tical Software, 40(1):1–29.
Wickham, H. (2017). tidyverse: Easily Install and Load the ’Tidyverse’. R package version 1.2.1.