R/plot.stslist.statd.R

In TraMineR: Trajectory Miner: a Sequence Analysis Toolkit

## ==============================
## PLot of the state distribution
## ==============================

plot.stslist.statd <- function(x, type = "d", cpal = NULL, ylab = NULL,
  yaxis = TRUE, xaxis = TRUE, xtlab = NULL, xtstep = NULL, tick.last = NULL,
    cex.axis = 1, space = 0, xlab = NULL, lwd=3.5, col="blue", ylim=NULL, cex.plot, ...) {

  TraMineR.check.depr.args(alist(cex.axis = cex.plot))

  if (!(type %in% c("d", "Ht", "dH"))){
    msg.stop("type can only be 'd', 'Ht', or 'dH'")
  }

  sep.ylab <- (isFALSE(yaxis) && (is.null(ylab) || !is.na(ylab)))
  cex.lab <- par("cex.lab")
  if ("cex.lab" %in% names(list(...))) cex.lab <- list(...)[["cex.lab"]]

	n <- attr(x,"nbseq")
	weighted <- attr(x, "weighted")
	if (weighted) {wlab <- "weighted "}
	else {wlab <- NULL}

	if (is.null(xtlab)) {xtlab <- attr(x,"xtlab")}

	if (is.null(xtstep)) {
		if (!is.null(attr(x,"xtstep"))) {xtstep <- attr(x,"xtstep")}
		## For sequence objects created with previous versions
		else {xtstep <- 1}
	}
	if(is.null(tick.last)){
		tick.last <- ifelse(!is.null(attr(x, "tick.last")), attr(x, "tick.last"), FALSE)
	}

	seql <- length(xtlab)

    ## for x axis
	x.lab.pos <- NULL
	tpos <- seq(1,seql, xtstep)
    if (tick.last & tpos[length(tpos)] < seql) tpos <- c(tpos,seql)

	for (p in tpos) {
		x.lab.pos <- c(x.lab.pos, (p-1)+((p-1)*space)+(0.5+space))
	}

    ## y label
    if(is.null(ylab)){
        if (type=="d")
			ylab <- paste("Rel. Freq. (",wlab,"n=",round(n,2),")",sep="")
        else
            ylab <- paste("Entropy (",wlab,"n=",round(n,2),")",sep="")
    }

    if (sep.ylab) {
        sylab <- ylab
        ylab <- NA
    }

    ## y limit for entropy
    y <- x$Entropy
    if (is.null(ylim)) {
        if (isTRUE(attr(x,"norm")))
            ylim <- c(0.025,1)
        else
            ylim <- c(0.025*1.1*max(y),1.1*max(y))
    }


	## State distribution plot
	if (type %in% c("d","dH")) {
		if (is.null(cpal))
			cpal <- attr(x,"cpal")

		#if (is.null(ylab))
		#	ylab <- paste("Freq. (",wlab,"n=",round(n,2),")",sep="")

		bp <- barplot(x$Frequencies,
			space=space,
			axes=FALSE,
			axisnames=FALSE,
			## cex.axis=cex.axis,
			## cex.names=cex.axis,
			col=cpal,
			ylab=ylab,
            xlab=xlab,
			...)

		## Plotting the x axis
 ##		x.lab.pos <- NULL
 ##		tpos <- seq(1,seql, xtstep)
 ##        if (tick.last & tpos[length(tpos)] < seql) tpos <- c(tpos,seql)
 ##
 ##		for (p in tpos) {
 ##			x.lab.pos <- c(x.lab.pos, (p-1)+((p-1)*space)+(0.5+space))
 ##		}

		if (xaxis) {
			axis(1, at=x.lab.pos, labels=xtlab[tpos], pos=-0.02, cex.axis=cex.axis)
		}

        if (type == "dH") {
            #if (is.null(ylab.r))
            #    ylab.r <- "Entropy index"
            par(new=TRUE)
            c1 <- 1 - (1-.5*space)/length(y)
            c2 <- .5
            bp <- as.vector(bp)[1:length(y)]
            #plot(x=as.vector(bp)[1:length(y)], y=y, type="n", axes=FALSE,
            plot(x=x.lab.pos, y=y, type="n", axes=FALSE,
                    ylim=ylim, xlab=NA, ylab=NA)
            #lines(x=c1*bp + c2, y=y, col=col, lwd=lwd)
            lines(x=c1*x.lab.pos + c2, y=y, col=col, lwd=lwd)
            #axis(4,at=seq(0,max(c(1,max(y))),.2))
            #mtext(ylab.r, line=3, side=4, cex = par("cex"))

        }
	}
	## Entropy index plot
	else if (type=="Ht") {
        #if (is.null(ylab.r))
        #    ylab.r <- paste("Entropy index (",wlab,"n=",round(n,2),")",sep="")

		plot(y,
			col=col,
			## frame.plot=TRUE,
			type="l",
			lwd=lwd,
			lty="solid",
			axes=FALSE,
			ylim=ylim,
			ylab=ylab,
			xlab=xlab,
			...)

		## Plotting the x axis
		if (xaxis) {
			#tpos <- seq(1,seql, xtstep)
			axis(1, at=tpos, labels=xtlab[tpos], pos=-0.02, cex.axis=cex.axis)
		}
	}

	##
    if (sep.ylab)
        title(ylab=sylab, line=1, cex.lab=cex.lab)
	if (is.null(yaxis) || yaxis)
		axis(2, cex.axis=cex.axis)

  if (type == 'Ht') return(invisible(x$Entropy))
  else if (type == 'd') return(invisible(x$Frequencies))
  else return(invisible(x))
}