
R version 4.3.3 (2024-02-29) -- "Angel Food Cake"
Copyright (C) 2024 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

  Natural language support but running in an English locale

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> library(ouch)
> 
> treedat <- data.frame(node=1:7,anc=NA,time=NA)
> treedat$anc[-1] <- c(1,2,2,1,5,5)
> treedat$time <- c(0,1,3,3,2,3,3)
> treedat$reg <- as.factor(c("a","a","a","b","b","a","b"))
> 
> treedat$data1[c(3,4,6,7)] <- c(-0.00133741562273541, -1.01339260146802, -0.0022452348823774, -1.00140366716313)
> treedat$data2[c(3,4,6,7)] <- c(0.0161586802663828, 1.01327463666717, 0.0125379175716403, 0.991540176709494)
> treedat$data3[c(3,4,6,7)] <- c(-4,-1,0,-2)
> treedat$data4[c(3,4,6,7)] <- c(0.5,1,-3,1)
> treedat$data5[c(3,4,6,7)] <- c(2,0,1,-2)
> 
> tree <- with(treedat,ouchtree(node,anc,time))
> print(tree)
  nodes ancestors times labels
1     1      <NA>     0      1
2     2         1     1      2
3     3         2     3      3
4     4         2     3      4
5     5         1     2      5
6     6         5     3      6
7     7         5     3      7
> show(tree)
  nodes ancestors times labels
1     1      <NA>     0      1
2     2         1     1      2
3     3         2     3      3
4     4         2     3      4
5     5         1     2      5
6     6         5     3      6
7     7         5     3      7
> x <- treedat[c(6,4,5,7,3,2,1),c("data1","data2")]
> btree <- brown(data=x,tree)
> print(btree)

call:
brown(data = x, tree = tree)
  nodes ancestors times labels        data1      data2
1     1      <NA>     0      1           NA         NA
2     2         1     1      2           NA         NA
3     3         2     3      3 -0.001337416 0.01615868
4     4         2     3      4 -1.013392601 1.01327464
5     5         1     2      5           NA         NA
6     6         5     3      6 -0.002245235 0.01253792
7     7         5     3      7 -1.001403667 0.99154018

sigma squared:
           [,1]       [,2]
[1,]  0.1888074 -0.1853472
[2,] -0.1853472  0.1819546

theta:
NULL
   loglik  deviance       aic     aic.c       sic       dof 
 13.07838 -26.15675 -16.15675  13.84325 -15.75955   5.00000 
> show(btree)

call:
brown(data = x, tree = tree)
  nodes ancestors times labels        data1      data2
1     1      <NA>     0      1           NA         NA
2     2         1     1      2           NA         NA
3     3         2     3      3 -0.001337416 0.01615868
4     4         2     3      4 -1.013392601 1.01327464
5     5         1     2      5           NA         NA
6     6         5     3      6 -0.002245235 0.01253792
7     7         5     3      7 -1.001403667 0.99154018

sigma squared:
           [,1]       [,2]
[1,]  0.1888074 -0.1853472
[2,] -0.1853472  0.1819546

theta:
NULL
   loglik  deviance       aic     aic.c       sic       dof 
 13.07838 -26.15675 -16.15675  13.84325 -15.75955   5.00000 
> htree <- hansen(data=x,tree,regimes=treedat["reg"],sqrt.alpha=c(1,0,1),sigma=c(0.1,0,0.1),fit=F)
> htree <- update(htree,sqrt.alpha=c(1,-0.1,1),sigma=c(0.5,0.5,1),fit=F)
> print(htree)

call:
hansen(data = data, tree = object, regimes = regimes, sqrt.alpha = sqrt.alpha, 
    sigma = sigma, fit = ..1)
  nodes ancestors times labels reg reg.1        data1      data2
1     1      <NA>     0      1   a     a           NA         NA
2     2         1     1      2   a     a           NA         NA
3     3         2     3      3   a     a -0.001337416 0.01615868
4     4         2     3      4   b     b -1.013392601 1.01327464
5     5         1     2      5   b     b           NA         NA
6     6         5     3      6   a     a -0.002245235 0.01253792
7     7         5     3      7   b     b -1.001403667 0.99154018

alpha:
     [,1]  [,2]
[1,]  1.0 -0.10
[2,] -0.1  1.01

sigma squared:
     [,1] [,2]
[1,] 0.25 0.25
[2,] 0.25 1.25

theta:
$data1
         a          b 
 0.1305478 -1.0640968 

$data2
         a          b 
-0.1571595  1.0860827 

    loglik   deviance        aic      aic.c        sic        dof 
 -2.433922   4.867845  24.867845 -48.465489  25.662260  10.000000 
> show(htree)

call:
hansen(data = data, tree = object, regimes = regimes, sqrt.alpha = sqrt.alpha, 
    sigma = sigma, fit = ..1)
  nodes ancestors times labels reg reg.1        data1      data2
1     1      <NA>     0      1   a     a           NA         NA
2     2         1     1      2   a     a           NA         NA
3     3         2     3      3   a     a -0.001337416 0.01615868
4     4         2     3      4   b     b -1.013392601 1.01327464
5     5         1     2      5   b     b           NA         NA
6     6         5     3      6   a     a -0.002245235 0.01253792
7     7         5     3      7   b     b -1.001403667 0.99154018

alpha:
     [,1]  [,2]
[1,]  1.0 -0.10
[2,] -0.1  1.01

sigma squared:
     [,1] [,2]
[1,] 0.25 0.25
[2,] 0.25 1.25

theta:
$data1
         a          b 
 0.1305478 -1.0640968 

$data2
         a          b 
-0.1571595  1.0860827 

    loglik   deviance        aic      aic.c        sic        dof 
 -2.433922   4.867845  24.867845 -48.465489  25.662260  10.000000 
> 
> htree <- hansen(
+   tree=tree,
+   data=treedat[c('data3','data4','data5')],
+   regimes=treedat['reg'],
+   sqrt.alpha=c(1,0.5,0.5,1,0,3),
+   sigma=c(0.5,-0.1,-0.25,1,0.5,1),
+   fit=F
+ )
> print(htree)

call:
hansen(data = treedat[c("data3", "data4", "data5")], tree = tree, 
    regimes = treedat["reg"], sqrt.alpha = c(1, 0.5, 0.5, 1, 
        0, 3), sigma = c(0.5, -0.1, -0.25, 1, 0.5, 1), fit = F)
  nodes ancestors times labels reg reg.1 reg.2 data3 data4 data5
1     1      <NA>     0      1   a     a     a    NA    NA    NA
2     2         1     1      2   a     a     a    NA    NA    NA
3     3         2     3      3   a     a     a    -4   0.5     2
4     4         2     3      4   b     b     b    -1   1.0     0
5     5         1     2      5   b     b     b    NA    NA    NA
6     6         5     3      6   a     a     a     0  -3.0     1
7     7         5     3      7   b     b     b    -2   1.0    -2

alpha:
     [,1] [,2] [,3]
[1,]  1.0 0.50 0.50
[2,]  0.5 1.25 0.25
[3,]  0.5 0.25 9.25

sigma squared:
       [,1]   [,2]    [,3]
[1,]  0.250 -0.050 -0.1250
[2,] -0.050  1.010  0.5250
[3,] -0.125  0.525  1.3125

theta:
$data3
        a         b 
-2.415739 -1.281838 

$data4
         a          b 
-0.9444604  0.6716862 

$data5
        a         b 
 1.597923 -1.084672 

   loglik  deviance       aic     aic.c       sic       dof 
-41.42685  82.85369 118.85369  21.13941 127.58201  18.00000 
> 
> htree <- hansen(
+   tree=tree,
+   data=treedat[c('data3','data4','data5')],
+   regimes=treedat['reg'],
+   sqrt.alpha=c(1,0.5,0.5,1,0,3),
+   sigma=c(0.5,-0.1,-0.25,1,0.5,1),
+   fit=F
+ )
> 
