Number of particles in clusters follow-up

Thanks again for the suggestion to use spatstat in Bio7. It looks like this is the way to go. I read a lot about Bio7 and spatstat in the last days and installed Bio7 (version 2.4) and Rserve (version 3.3.1) and spatstat in the library today. Everything works how it should. I can open my image in imagej export the image size and matrix but when I try to use plot() I get this message

 *** caught segfault ***
address 0x110, cause 'memory not mapped'

Traceback:
 1: tiff(filename, width = 480, height = 480)
 2: .bio7Device()
 3: par()
 4: resolve.defaults(list(...), par())
 5: plot.owin(BB, type = "n", add = add, main = blankmain, show.all = show.all,     cex.main = cex.main)
 6: plot(BB, type = "n", add = add, main = blankmain, show.all = show.all,     cex.main = cex.main)
 7: plot.ppp(X)
 8: plot(X)
 9: doTryCatch(return(expr), name, parentenv, handler)
10: tryCatchOne(expr, names, parentenv, handlers[[1L]])
11: tryCatchList(expr, classes, parentenv, handlers)
12: tryCatch(expr, error = function(e) {    call <- conditionCall(e)    if (!is.null(call)) {        if (identical(call[[1L]], quote(doTryCatch)))             call <- sys.call(-4L)        dcall <- deparse(call)[1L]        prefix <- paste("Error in", dcall, ": ")        LONG <- 75L        msg <- conditionMessage(e)        sm <- strsplit(msg, "\n")[[1L]]        w <- 14L + nchar(dcall, type = "w") + nchar(sm[1L], type = "w")        if (is.na(w))             w <- 14L + nchar(dcall, type = "b") + nchar(sm[1L],                 type = "b")        if (w > LONG)             prefix <- paste0(prefix, "\n  ")    }    else prefix <- "Error : "    msg <- paste0(prefix, conditionMessage(e), "\n")    .Internal(seterrmessage(msg[1L]))    if (!silent && identical(getOption("show.error.messages"),         TRUE)) {        cat(msg, file = stderr())        .Internal(printDeferredWarnings())    }    invisible(structure(msg, class = "try-error", condition = e))})
13: try((plot(X)))
14: doTryCatch(return(expr), name, parentenv, handler)
15: tryCatchOne(expr, names, parentenv, handlers[[1L]])
16: tryCatchList(expr, classes, parentenv, handlers)
17: tryCatch(expr, error = function(e) {    call <- conditionCall(e)    if (!is.null(call)) {        if (identical(call[[1L]], quote(doTryCatch)))             call <- sys.call(-4L)        dcall <- deparse(call)[1L]        prefix <- paste("Error in", dcall, ": ")        LONG <- 75L        msg <- conditionMessage(e)        sm <- strsplit(msg, "\n")[[1L]]        w <- 14L + nchar(dcall, type = "w") + nchar(sm[1L], type = "w")        if (is.na(w))             w <- 14L + nchar(dcall, type = "b") + nchar(sm[1L],                 type = "b")        if (w > LONG)             prefix <- paste0(prefix, "\n  ")    }    else prefix <- "Error : "    msg <- paste0(prefix, conditionMessage(e), "\n")    .Internal(seterrmessage(msg[1L]))    if (!silent && identical(getOption("show.error.messages"),         TRUE)) {        cat(msg, file = stderr())        .Internal(printDeferredWarnings())    }    invisible(structure(msg, class = "try-error", condition = e))})
18: try(try((plot(X))))
19: print(try(try((plot(X)))))
20: eval(expr, envir, enclos)
21: eval(expr, pf)
22: withVisible(eval(expr, pf))
23: evalVis(expr)
24: capture.output(print(try(try((plot(X))))))
25: paste(capture.output(print(try(try((plot(X)))))), collapse = "\n")

Possible actions:
1: abort (with core dump, if enabled)
2: normal R exit
3: exit R without saving workspace
4: exit R saving workspace
Selection: 

I am not sure if I forgot to install another package but I can’t find any solution for this. When I run only R (without Bio7) I can use the command “plot()”.

Maybe you can help me again,
Joko

1 Like

Hello that’s strange. Which Operating System do you use?

This seems to be a plot device error.

If you use macosx try to open a device with the command quartz() before you plot.

Please note that you can change to the native connection and native r device without loosing the current workspace by hitting the Start Rserve Action in Bio7 again. This works forth and back and allows you to use R from within Bio7 without Rserve but with the same workspace.

See:

What happens if you plot in the native connection?

Momentarily I’m on a conference until thursday without a Laptop to reproduce the error if you name me the OS.

1 Like

Thanks I’ll try it tomorrow morning when I’m back in the office.
I use the latest MacOS.

Í recommend to install the R Mac GUI. It comes with all Devices installed.

https://cran.r-project.org/bin/macosx/

I hope you also installed the special Rserve package as described in the manual. It has to run in cooperative mode for a shared R workspace.

Yes, that is the R Mac GUI I’m running.

I think Rserve is the problem. When I run plot() in the native connection everything works quite well.

I had some problems to install Rserve because I got this error:

Error in install.packages("/Users/_username_/Downloads/Rserve_1.8-4_Mac_cooperative.tar",  : 
  type == "both" cannot be used with 'repos = NULL'

therefore I unpacked the file in the library directory and R shows me Rserve as installed package in the “R package installer”

When I start R in Bio7 I get this and assume that Rserve is running:

R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-apple-darwin13.4.0 (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.

> options(max.print=5000)
> cat("
+ ")

> options(device='cairo')
> library(Rserve)
> run.Rserve()
-- running Rserve in this R session (pid=12723), 1 server(s) --
(This session will block until Rserve is shut down)
Rserve connection established !
1 Like

Please download and install Rserve again. It is important that you don’t unzip the file which Safari does.

Here is how you can avoid the uncompressing.

https://discussions.apple.com/thread/1483114?start=0

You have to install the tgz, See the exact command:

http://bio7.org/manual/Main.html#toc-Subsubsection-2.1.2.3

Then the install command will work and you get a valid Installation.

Before remove your unziped Rserve!!

If that doesn’t work what happens if you change the plot Devices?

See:
http://bio7.org/?p=2112

1 Like

Thank you. I installed everything again and now Rserve is working how it supposed to work.
I start to go through the first examples and tried to run your scripts on my images.
I have some more questions regarding this:

  1. When I transfer my particle positions to R and plot them the position is mirrored because ImageJ sets the origin of xy in the top left corner. I think therefore you use plot(W, ylim=rev(W$yrange),main=NULL). That works for plot(X, add=T) but the density plot is still displaced.

  2. This was just a small ROI of my whole Image. Is there a limit how many particles or image sizes Bio7 can handle?

  3. Is it possible get a better resolution if it comes to more particles like in this case with 8004 particles?


  4. Like you see in my images I have just a small ROI with high particle density. I assume if I perform the L-function it will always calculate the distribution of the whole image and not just my selection and will always indicate some sort of clustering. Is it possible to import something like a selection from ImageJ and use just this are for the L-function?

I am sorry for all the questions, but spatstat is completely new to me.

Thanks for all your help and have a nice day,
Joko

1 Like
  1. Normally if you add the density plot like my posted examples add=T then the plot should be flipped. However I
    will reproduce it on a Mac when I’m back.

  2. The limits of points spatstat and therefore Bio7 can handle are described here:

https://spatstat.github.io/FAQ.html

If you have more RAM you can increase the size for Java and Bio7, See:

http://bio7.org/?page_id=37

R is independent of this settings and uses all available memory. So a sensitive Balance between Java and R if necessary can be applied if the memory use is at the limit.

The limit for Images is the adress space of Java using integers for indexing, see:

http://imagej.1557.x6.nabble.com/what-s-the-maximum-image-size-ImageJ-can-open-and-process-td3685418.html

However you could also workaround this by sending chunks to R if that would be necessary.

3 You can play with the parameter adjust to define the density Kernel in the density plot. You can change the plot size in the Bio7 plot preferences.

4 In spatstat the window is the defining Area and you can send selection coordinates from Bio7 to R for the use of defining the window. You could also do this in spatstat when you define your window. I can send you links later but please consult the spatstat documentation.

1 Like

Here you find free chapters of the spatstat book with explanations: See chapter 3:

https://spatstat.github.io/book/

I can reproduce the flipped coordinates on MacOSX. Very strange. I think this is a spatstat bug.

If you use select the pdf device (Context menu-> Plot Preferences - Acrobat reader installed!) then everything is displayed correctly.

I opened a bug here:

2 Likes

Thank you for all the help and also the spatstat bug.

I am sorry for the late reply, this is just a small project I am working on at the moment. Finally, I wrote some code for my cluster analysis. My current issue is, that I have 48 images to analyse and need to get there ppp (adjust-Threshold; Analyze Particle) to R and then run my R script to plot some of the patterns and save the results of my neighbor analysis in a csv file. Later I need to do a sensitivity analysis and want to avoid opening ImageJ every time.

library(spatstat)
X<- ppp(Particles$X, Particles$Y, c(0,imageSizeX), c(0,imageSizeY))
nnneighbor <- applynbd(X, R = 30, function(Y, ...){npoints(Y)}, exclude = TRUE)
result<- cbind(Particles, nnneighbor)
write.csv(result,"/Users/xxx/Documents/Cluster/Cluster.csv")

#show histogram and save it
hist(result$nnneighbor, breaks=100)

jpeg('/Users/xxx/Documents/Cluster/Cluster_R=30.jpg')
hist(result$nnneighbor, breaks=100)
dev.off()

P<- subset(result, nnneighbor>5)

jpeg("/Users/xxx/Documents/Cluster/27-094-L.jpg")
plot(result$X, result$Y,  xlim=c(0,30000), ylim=c(0,30000))  #plot all ppp
title(main='all ppp',xlab='x',ylab='y')
dev.off()

jpeg("/Users/xxx/Documents/Cluster/27-094-L_neighbors.jpg")
plot(P$X,P$Y, xlim=c(0,30000), ylim=c(0,30000)) #plot only nnneighbor>5
title(main='>5 neighbors',xlab='x',ylab='y')
dev.off()

How do I create an easy batch mode process in Bio7 to get my ppp information to R? I just found the ClusterMultipleFiles.bsh example and can’t make sense out of the code because I don’t know Bash.
Is there a simple example for this or something like a macro recorder (like in ImageJ) to create a macro?

1 Like

There are some ways to automate this in Java, Jython, BeanShell and Groovy scripts whatever one prefers.

And then you can methodically for instance import an image sequence as a virtual stack (disk resident) and then loop and analyze all stack slices (images). See the links at the end.

However I decided to create a simple BeanShell script (it is not a shell script) because ImageJ supports this scripting type (select it in the recorder) in the ImageJ macro recorder to record all ImageJ actions for a script.

BeanShell is a very old but easy to use Java scripting language still popular because you can also execute Java source up to version 1.5.

Just create this script in the project folder along with your R script. You R script will be called at the end of this script.

In the file dialog select all images you want to analyze (Ctrl+Mouse-Click or Shift+Mouse-Click).

Note, that all measured particle values are stored in seperate dataframes which you can then edit.

Change the script to your needs (your R script only save the last image values at the moment. You can e.g., iterate over all dataframes, etc.). For a faster speed I do not display the images.

it is also possible to put all R comands in the BeanShell script if you like. But I created a special method in the Bio7 API RServe.callRScriptSilent to call a R script from within another Java scripting language to avoid a quoting hell of commands. I possible collect the values first (as I do below) and then plot and save at the end (change the script!) so you don’t have to call the R script unnecessary in a loop. This will make a hughe speed difference.

/*A BeanShell script to transfer a Particle Analysis to R!*/

import ij.IJ;
import ij.ImagePlus;
import ij.process.ImageProcessor;
import org.rosuda.REngine.Rserve.RConnection;
import org.rosuda.REngine.Rserve.RserveException;
import com.eco.bio7.image.r.IJTranserResultsTable;
import com.eco.bio7.rbridge.RServe;
import com.eco.bio7.rbridge.RState;

if (RServe.isAliveDialog()) {
    if (RState.isBusy() == false) {
        String[] files = Bio7Dialog.openMultipleFiles();//Open files from multiple selections(with Shift+Mouse-Click or Ctrl+Mouse-Click).
        if (files != null) {
            RConnection c = RServe.getConnection();
            for (int i = 0; i < files.length; i++) {
                System.out.println(Bio7Dialog.getCurrentPath() + "/" + files[i]);
                ImagePlus imp = IJ.openImage(Bio7Dialog.getCurrentPath() + "/" + files[i]);
               /*Open  image data with the ImageJ API if needed!*/
                //imp.show();
                /*Here we execute some BeanShell commands which can be produced by the ImageJ macro recorder!*/
                IJ.run(imp, "8-bit", "");
                IJ.setAutoThreshold(imp, "Default");
                IJ.run(imp, "Analyze Particles...", "display clear");
                
                /*Transfer the results table! We add an empty second argument (ImageJ macro) to simply transfer the results table. The default API makes a particle analysis!*/
                IJTranserResultsTable.runParticleAnalysis(c, "");

                /*Here we create new individual dataframe with the loop count as prefix for the name to create several measurements!*/
                c.eval("assign(paste0('Particles'," + i + ",sep=''),Particles)");

                /*We need to assign the image size. They are not automatically transferred!*/
                c.eval("imageSizeX<-" + imp.getWidth() + "");
                c.eval("imageSizeY<-" + imp.getHeight() + "");

                /*Close the image to free up space if needed!*/
               // imp.close();

                /*Call the R script in the project (relative path!)*/
                RServe.callRScriptSilent("/LoopImagesR/loop.R");

            }
        }
    }
}

Some related scripts can be found here:

Iterate over a stack example:

https://github.com/Bio7/bio7_examples/blob/master/Bio7_Script_Examples/Cluster_Images_Examples/ClusterRGBStack.bsh

Open images and transfer the values to the Bio7 spredsheet (store values in Excel using the spreadsheat actions)

https://github.com/Bio7/bio7_examples/blob/master/Bio7_Script_Examples/Cluster_Images_Examples/ClusterMultiple%20Files.bsh

Apropos, if you use Java for this (ImageJ commands also recordable) you could also create a nice interface for this kind of analysis by using the powerful SWT WindowBuilder of Eclipse or the embedded JavaFX SceneBuilder.

1 Like

Thanks again for the help. I tested the bsh. script and it works fine with two exceptions:

  1. I have a batch of 47 images and after 23 of them Bio7 isn’t responding any more and I have to force quit it. Can this be a memory issue? I have the same issue when I run the batch mode in imageJ with large files.
  2. Can I give the name of the currently open file from BeanShell to my R script to name my plots and csv file? At the moment I just renumber (from 1-47) because I don’t know how to get the title of files[i] from BeanShell to R.

Ok I solved issue 1). Image no 23 was too big. I cropped the png image from 1.1MB to 0.85MB and now it works

Issue 1: You could also increase the available RAM by increasing the memory for Bio7, from the Bio7 FAQ:

How can i increase the available RAM for Bio7?

Windows and Linux:

For an increased Java heap space open the Bio7.ini file in the install
directory of Bio7. In the file you can change the default memory
settings e.g. the initial heap size -Xms and the maximum heap space
-Xmx.

MacOSX:

For an increased Java heap space open the Bio7 package (context menu if
you click on the icon) then go to Contents->MacOS and open the
Bio7.ini file with a texteditor. In the file you can change the default
memory settings e.g. the initial heap size -Xms and the maximum heap
space -Xmx

I would suggest you spent half of the RAM of the OS for Java to leave some memory for R. You can monitor the increase of used RAM quite easily with your OS systems monitor to get a feeling for the memory consumption of the Java and R parts, see:

An intelligent solution is using both initial heap Xms and maximum heap Xmx. Only if the memory is required Java uses the maximum heap (easy to monitor).

You could also close the images in the loop to free space with the command imp.close(); instantly. But the variable is reused in the loop and no further big Java object is created in the loop.

A hughe speed difference is to save all results outside the loop at the end: As I already mentioned all image results are dynamically stored in the R workspace which you can save at the end in one shot. Just comment out the command RServe.callRScriptSilent("/LoopImagesR/loop.R");
to see the difference. I also think that a looped R save operation leads to crashes on the R side because of locked files on the R files. But this has to be tested more.

Issue 2: You could easily transfer the image name to the R workspace, e.g., with:

c.assign("currentImageName",imp.getShortTitle());
c.assign("currentImageNameLong", imp.getTitle() );

See: http://javadoc.imagej.net/ImageJ1/ij/ImagePlus.html

Or create all measurements named in the workspace within the loop with:

/*Here we create new individual dataframe with the loop count and image name as prefix for the name to create several measurements!*/   
c.assign("tempTitle","_"+imp.getShortTitle()+"_"+i); //transfer a Java string to the R workspace!
c.eval("assign(paste0('Particles',tempTitle,sep=''),Particles)"); //Concatenate to whatever title!

In the last line I uses the ‘assign’ command on the R side to concatenate a R string ‘Particles’ with the transferred Java String.

1 Like

By the way on StackOverflow I answered the problem how to make you massive amount of points more visible:

http://stackoverflow.com/questions/39258578/how-do-i-improve-the-resolution-of-the-plotdensity-in-spatstat

Thanks again, c.assaign(“currentImageName”… really helped a lot and I will change the memory settings soon after I imported all of my ImageJ plugins because I’ll try to do all my future ImageJ work in Bio7 now.
I’ll let you know if I have still memory issues while using ImageJ batch mode in Bio7 and will perform a sensitivity analysis of my data at the weekend if I finish my code tomorrow.

And of course, you can make this topic public again.

Best,
Joko

I finished my cluster analysis and just need to transfer my particles in the clusters (e.g. defined as particles with >5 neighbors in a 30px distance) as a selection or mask back into my original images to have a visual control for review and presentation.
My current workflow looks like this:

  1. a stack of images (46 with ~1 GB each)
    → thresholding and particle analysis
    → create mask
    → save mask as png
dir1 = getDirectory("Choose Source Directory ");
dir2 = getDirectory("Choose Source Directory ");
list = getFileList(dir1);
setBatchMode(true);
for (n=0; n<list.length; n++) {
    showProgress(n+1, list.length);
    open(dir1+list[n]);

run("Clear Outside");
run("Threshold Colour");
// Threshold Colour v1.12a------
// Autogenerated macro, single images only!
// G. Landini 27/Sep/2010.
//
min=newArray(3);
max=newArray(3);
filter=newArray(3);
a=getTitle();
run("HSB Stack");
run("Convert Stack to Images");
selectWindow("Hue");
rename("0");
selectWindow("Saturation");
rename("1");
selectWindow("Brightness");
rename("2");
min[0]=0;
max[0]=255;
filter[0]="pass";
min[1]=45;
max[1]=255;
filter[1]="pass";
min[2]=0;
max[2]=255;
filter[2]="pass";
for (i=0;i<3;i++){
  selectWindow(""+i);
  setThreshold(min[i], max[i]);
  run("Convert to Mask");
  if (filter[i]=="stop")  run("Invert");
}
imageCalculator("AND create", "0","1");
imageCalculator("AND create", "Result of 0","2");
for (i=0;i<3;i++){
  selectWindow(""+i);
  close();
}
selectWindow("Result of 0");
close();
selectWindow("Result of Result of 0");
rename(a);
// Threshold Colour ------------
run("8-bit");
setAutoThreshold("Default");
//run("Threshold...");
//setThreshold(0, 254);
setOption("BlackBackground", false);
run("Convert to Mask");
setThreshold(1, 255);
run("Analyze Particles...", "size=10-1400 show=Masks display exclude summarize add");
run("Watershed");
run("Analyze Particles...", "size=10-400 circularity=0.40-1.00 show=Masks display exclude summarize add");

saveAs("png", dir2+list[n]);
run("Close All");
}
  1. load png in Bio7
    → create ppp in spatstat
/*A BeanShell script to transfer a Particle Analysis to R!*/

import ij.IJ;
import ij.ImagePlus;
import ij.process.ImageProcessor;
import org.rosuda.REngine.Rserve.RConnection;
import org.rosuda.REngine.Rserve.RserveException;
import com.eco.bio7.image.r.IJTranserResultsTable;
import com.eco.bio7.rbridge.RServe;
import com.eco.bio7.rbridge.RState;


if (RServe.isAliveDialog()) {
    if (RState.isBusy() == false) {
        String[] files = Bio7Dialog.openMultipleFiles();//Open files from multiple selections(with Shift+Mouse-Click or Ctrl+Mouse-Click).
        if (files != null) {
            RConnection c = RServe.getConnection();
            for (int i = 0; i < files.length; i++) {
                System.out.println(Bio7Dialog.getCurrentPath() + "/" + files[i]);
                ImagePlus imp = IJ.openImage(Bio7Dialog.getCurrentPath() + "/" + files[i]);
               /*Open  image data with the ImageJ API if needed!*/
                //imp.show();
                /*Here we execute some BeanShell commands which can be produced by the ImageJ macro recorder!*/
                IJ.run(imp, "8-bit", "");
                IJ.setAutoThreshold(imp, "Default");
                IJ.run(imp, "Analyze Particles...", "display clear");
                
                /*Transfer the results table! We add an empty second argument (ImageJ macro) to simply transfer the results table. The default API makes a particle analysis!*/
                IJTranserResultsTable.runParticleAnalysis(c, "");

                /*Here we create new individual dataframe with the loop count as prefix for the name to create several measurements!*/
                c.eval("assign(paste0('Particles'," + i + ",sep=''),Particles)");

                /*We need to assign the image size. They are not automatically transferred!*/
                c.eval("imageSizeX<-" + imp.getWidth() + "");
                c.eval("imageSizeY<-" + imp.getHeight() + "");

                /*Close the image to free up space if needed!*/
               // imp.close();
                c.assign("currentImageName",imp.getShortTitle());
                
                /*Call the R script in the project (relative path!)*/
                RServe.callRScriptSilent("/Clustering/nneighbor.R");
        
            }
        }
    }
}

→ get list of particles in clusters and plot them for a first look if I selected the right region and particles

> library(spatstat)

> dirname <- "/Users/xxx/Documents/Cluster/neighbor"
> X<- ppp(Particles$X, Particles$Y, c(0,imageSizeX), c(0,imageSizeY))
> nnneighbor <- applynbd(X, R = 30, function(Y, ...){npoints(Y)}, exclude = TRUE)
> result<- cbind(Particles, nnneighbor)
> write.csv(result,file.path(dirname,"/","results","/",x, paste(currentImageName,".csv", sep = "")))
> }

> P<- subset(result, nnneighbor>6)

> jpeg(file.path(dirname,"/","nnneighbor>6", paste(currentImageName,".jpeg", sep = "")))
> plot(P$X,P$Y, xlim=c(0,30000), ylim=c(0,30000)) #plot only nnneighbor>6
> title(main='>6 neighbors',xlab='x',ylab='y')
> dev.off()

I didn’t integrate my ImageJ macro in Bio7 yet because I run out of memory every 2-3 images (in ImageJ and also in Bio7) no matter how I change the memory settings because the memory accumulates with every new image. But that’s a minor issue. I can also do this image by image (for the 46 images), although it would be nice to find a solution to process the whole batch because I will have a lot more images in the future (I addressed this issue in another topic in this forum).

I thought I can use the ROI manager to delete/deselect all of the particles in my first selection that have less than 7 neighbors because I get a list of IDs (I think the IDs at the end match the IDs in the ROI manager) with x/y position and the number of neighbors.
I am not sure how to implement this and if that is possible. Do you have any idea what the best approach for this is?

Have a nice weekend,
Joko

  1. Could you send me an example image so that I can analyze the memory consumption?
    Possibly there must be a memory leak on the ImageJ Java side or caused by an analysis step in your specific macro. Please note that the ImageJ1 API is used on the Bio7 side by default to open images so it might not be related to the opened bug report since you can reproduce it on FIJI and Bio7. Apropos you can make an easy test what happens if you just open and close all images in a short macro. Do you also run out of memory? If not I would do a step by step analysis (debugging) to find the memory leak.

  2. Do you mean the id of a R dataframe column from a spatstat analysis? It is definitely possible to select the results of an R analysis and map it to the ROI manager to select ROI’s with a certain ID. Recently I helped someone to create a bidrectional GUI to create an ImageJ->R->ImageJ analysis.