kuwisdelu / Cardinal

Mass spectrometry imaging toolbox
http://www.cardinalmsi.org
36 stars 14 forks source link

Missing TRUE/FALSE value in: if ( all(class==pred$class) ) break #11

Closed gtluu closed 4 years ago

gtluu commented 4 years ago

https://github.com/kuwisdelu/Cardinal/blob/f31dc7dc6e0cef3c40708cec261e3b54fd239424/R/stats2-spatialShrunkenCentroids.R#L215

Hi Kylie,

I seem to have an issue in in spatialShrunkenCentroids (Cardinal 2.4.0) where the above line produces an error when the result of the statement all(class==pred$class)) == NA. I believe this is occurring in my dataset when a pixel/pixels cannot be classified.

I have found a workaround in the following code for now by allowing comparisons to include vectors of NA values as a workaround that seems to produce the same results, and when a pixel cannot be classified, it is simply essentially removed from the segmentation map (appears as a white pixel in the image). This seems to work, though I'm curious as to whether or not this causes any erroneous results in SSC.

Thanks!

.spatialShrunkenCentroids2_cluster <- function(x, r, k, s, mean, class, weights,
                                                dist, iter.max, BPPARAM)
{
    iter <- 1
    init <- TRUE
    # suppress progress in inner parallel loop
    progress <- getOption("Cardinal.progress")
    options(Cardinal.progress=FALSE)
    on.exit(options(Cardinal.progress=progress))
    # suppress messages
    verbose <- getOption("Cardinal.verbose")
    # cluster the data
    while ( iter <= iter.max ) {
        options(Cardinal.verbose=FALSE) # suppress messages
        class <- factor(as.integer(droplevels(class)))
        if ( nlevels(class) == 1L )
            .stop("singular cluster configuration, try a different seed")
        fit <- .spatialShrunkenCentroids2_fit(x,
            s=s, mean=mean, class=class, BPPARAM=BPPARAM)
        pred <- .spatialShrunkenCentroids2_predict(x,
            r=r, class=class, weights=weights, centers=fit$centers,
            sd=fit$sd, priors=NULL, dist=dist, init=init, BPPARAM=BPPARAM)
        # MODIFIED CODE
    compareClasses <- function(c1, c2) {
        same <- all(c1==c2)
        if (!is.na(same) & same==TRUE) {
            return(TRUE)
        } else if (is.na(same) | same==FALSE) {
            return(FALSE)
        }
    }
        if ( compareClasses(c1=class, c2=pred$class) )
            break
        # END MODIFIED CODE
        class <- pred$class
        init <- pred$init
        iter <- iter + 1
        options(Cardinal.verbose=verbose) # print progress now
        .message(".", appendLF=FALSE)
    }
    options(Cardinal.verbose=verbose) # restore
    .message(" ")
    list(class=pred$class, probability=pred$probability, scores=pred$scores,
        centers=fit$centers, statistic=fit$statistic, sd=fit$sd)
}

.spatialShrunkenCentroids2_predict <- function(x, r, class, weights,
                            centers, sd, priors, dist, init, BPPARAM)
{
    # suppress progress in inner parallel loop
    progress <- getOption("Cardinal.progress")
    options(Cardinal.progress=FALSE)
    on.exit(options(Cardinal.progress=progress))
    # calculate spatial discriminant scores
    fun <- function(xbl) {
        .spatialScores(xbl, centers=centers, weights=attr(xbl, "params"),
            neighbors=attr(xbl, "neighbors"), sd=sd + s0)
    }
    s0 <- median(sd)
    scores <- spatialApply(x, .r=r, .fun=fun, .blocks=TRUE,
        .simplify=.rbind_and_reorder, .init=init, .dist=dist,
        .params=weights, BPPARAM=BPPARAM)
    if ( isTRUE(init) ) {
        init <- attr(scores, "init")
        attr(scores, "init") <- NULL
    }
    if ( is.null(priors) && !is.null(class) )
        priors <- table(class) / sum(table(class))
    scores <- scores - 2 * log(rep(priors, each=ncol(x)))
    colnames(scores) <- colnames(centers)
    # calculate posterior class probabilities
    pfun <- function(scores) {
        s <- exp(-0.5 * scores)
        pmax(s / rowSums(s), 0)
    }
    probability <- pfun(scores)
    colnames(probability) <- colnames(centers)
    # calculate and return new class assignments
    # MODIFIED LINE
    class <- factor(apply(probability, 1, function(x) ifelse(length(x), which.max(x), NA)),
        levels=seq_len(ncol(centers)), labels=colnames(centers))
    list(scores=scores, probability=probability, class=class, init=init)
}

EDIT: updated to correct modified code. also found at: https://github.com/gtluu/Cardinal/blob/7a6a768cbb741cfcaf864e97ab6e0a063c093e4a/R/stats2-spatialShrunkenCentroids.R#L219

kuwisdelu commented 4 years ago

This should now be fixed by 4cf89e0