Learn R Programming

PST (version 0.94.1)

pdist: Compute probabilistic divergence between two PST

Description

Compute probabilistic divergence between two PST

Usage

# S4 method for PSTf,PSTf
pdist(x,y, method="cp", l, ns=5000, symetric=FALSE, output="all")

Value

If ouput="all", a vector containing the divergence value for each generated sequence, if output="mean", the mean, i.e. expected value which is the divergence between models.

Arguments

x

a probabilistic suffix tree, i.e., an object of class "PSTf" as returned by the pstree, prune or tune function.

y

a probabilistic suffix tree, i.e., an object of class "PSTf" as returned by the pstree, prune or tune function.

method

character. Method for computing distances. So far only one method is available.

l

integer. Length of the sequence(s) to generate.

ns

integer. Number sequences to generate.

symetric

logical. If TRUE, the symetric version of the measure is returned, see details.

output

character. See value.

Author

Alexis gabadinho

Details

The function computes a probabilistic divergence measure between PST \(S_{A}\) and \(S_{B}\) based on the measure originally proposed in Juang-1985 and Rabiner-1989 for the comparison of two (hidden) Markov models \(S_{A}\) and \(S_{B}\) $$ d(S_{A}, S_{B})=\frac{1}{\ell} [\log P^{S_{A}}(x)-\log P^{S_{B}}(x)]=\frac{1}{\ell}\log \frac{P^{S_{A}}(x)}{P^{S_{B}}(x)} $$ where \(x=x_{1}, \ldots, x_{\ell}\) is a sequence generated by model \(S_{A}\), \(P^{S_{A}}(x)\) is the probability of \(x\) given model \(S_{A}\) and \(P^{S_{B}}(x)\) is the probability of \(x\) given model \(S_{B}\). The ratio between the two sequence likelihoods measures how many times the sequence \(x\) is more likely to have been generated by \(S_{A}\) than by \(S_{2}\).

As the number \(n\) of generated sequences on which the measure is computed (or the length of a single sequence) approaches infinity, the expected value of \(d(S_{A}, S_{B})\) converges to \(d_{KL}(S_{A}, S_{B})\) Falkhausen-1995, He-2000, the Kullback-Leibler (KL) divergence (also called information gain) used in information theory to measure the difference between two probability distributions.

The pdist function uses the following procedure to compute the divergence between two PST:

  • generate a ransom sample of \(n\) sequences (of length \(\ell\)) with model \(S_{A}\) using the generate method

  • predict the sequences with \(S_{A}\) and with \(S_{B}\)

  • compute $$ d_{i}(S_{A}, S_{B})=\frac{1}{\ell} [\log P^{S_{A}}(x_{i})-\log P^{S_{B}}(x_{i}))], \; i=1, \ldots, n $$

  • the expected value $$ E(d(S_{A}, S_{B})) $$ is the divergence between models \(S_{A}\) and \(S_{B}\) and is estimated as $$ \hat{E}(d(S_{A}, S_{B}))=\frac{1}{n} \sum_{i=1}^{n} d_{i}(S_{A}, S_{B}) $$

For more details, see Gabadinho 2016.

References

Gabadinho, A. & Ritschard, G. (2016). Analyzing State Sequences with Probabilistic Suffix Trees: The PST R Package. Journal of Statistical Software, 72(3), pp. 1-39.

Juang, B. H. and Rabiner, L. R. (1985). A probabilistic distance measure for hidden Markov models. ATT Technical Journal, 64(2), pp. 391-408.

Rabiner, L. R. (1989). A tutorial on hidden Markov models and selected applications in speech recognition. Proceedings of the IEEE, 77(2), pp. 257-286.

Examples

Run this code
## activity calendar for year 2000
## from the Swiss Household Panel
## see ?actcal
data(actcal)

## selecting individuals aged 20 to 59
actcal <- actcal[actcal$age00>=20 & actcal$age00 <60,]

## defining a sequence object
actcal.lab <- c("> 37 hours", "19-36 hours", "1-18 hours", "no work")
actcal.seq <- seqdef(actcal,13:24,labels=actcal.lab)

## building a PST segmented by age group
gage10 <- cut(actcal$age00, c(20,30,40,50,60), right=FALSE,
	labels=c("20-29","30-39", "40-49", "50-59"))

actcal.pstg <- pstree(actcal.seq, nmin=2, ymin=0.001, group=gage10)

## pruning
C99 <- qchisq(0.99,4-1)/2
actcal.pstg.opt <- prune(actcal.pstg, gain="G2", C=C99)

## extracting PST for age group 20-39 and 30-39
g1.pst <- subtree(actcal.pstg.opt, group=1)
g2.pst <- subtree(actcal.pstg.opt, group=2)

## generating 5000 sequences with g1.pst 
## and computing 5000 distances
dist.g1_g2 <- pdist(g1.pst, g2.pst, l=11)
hist(dist.g1_g2)

## the probabilistic distance is the mean
## of the 5000 distances
mean(dist.g1_g2)

Run the code above in your browser using DataLab