Skip to content

Commit

Permalink
Calculate Mean, Min and Max for 5D Images
Browse files Browse the repository at this point in the history
Statistical tests that calculate the mean, min, and max for N5 files loaded into images.
  • Loading branch information
AvocadoMoon committed Jul 26, 2024
1 parent 6a0cf5f commit b20a80e
Show file tree
Hide file tree
Showing 2 changed files with 127 additions and 29 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -10,8 +10,8 @@ public class N5DataSetFile {
public final double[][] histMin;
public final double[][] histAverage;
public final double totalArea;
public final double testDomainArea;
public N5DataSetFile(String uri, String[] variables, HashMap<Integer, String> mask, double[][] histMax, double[][] histMin, double[][] histAverage, double totalArea, double testDomainArea){
public final double[] testDomainArea;
public N5DataSetFile(String uri, String[] variables, HashMap<Integer, String> mask, double[][] histMax, double[][] histMin, double[][] histAverage, double totalArea, double[] testDomainArea){
this.uri = uri;
this.histMax = histMax;
this.histMin = histMin;
Expand All @@ -23,7 +23,7 @@ public N5DataSetFile(String uri, String[] variables, HashMap<Integer, String> ma
}

public static N5DataSetFile[] alphaTestFiles(){
N5DataSetFile frapSimulationResultsMasked = new N5DataSetFile("https://vcell-dev.cam.uchc.edu/n5Data/ezequiel23/c607b779af9481f.n5?dataSetName=6262029569",
N5DataSetFile frapSimulationResultsMasked = new N5DataSetFile("https://vcell-dev.cam.uchc.edu/n5Data/ezequiel23/c607b779af9481f.n5?dataSetName=4457702594",
new String[]{"Dex"},
new HashMap<Integer, String>(){{put(1, "Cyt"); put(0, "Ec");}},
new double[][]{{10.0, 9.990392675155721, 9.83580092714469, 9.520539931524715, 9.162150060086567, 8.82335160436397, 8.523689113752786, 8.265381795870683, 8.044751960699015, 7.856809648125466,
Expand All @@ -35,9 +35,30 @@ public static N5DataSetFile[] alphaTestFiles(){
6.728509585652239, 6.7285095856523345, 6.728509585652324, 6.72850958565233, 6.728509585652335, 6.728509585652381, 6.728509585652381, 6.728509585652362, 6.728509585652365,
6.728509585652356, 6.7285095856523185, 6.728509585652278, 6.7285095856522865, 6.728509585652241, 6.728509585652113, 6.728509585652124, 6.7285095856521355, 6.728509585652127, 6.728509585652108}},
484,
314.64);
new double[]{0,314.64});

return new N5DataSetFile[]{frapSimulationResultsMasked};
N5DataSetFile anns5DTIRFSimulation = new N5DataSetFile("https://vcell-dev.cam.uchc.edu/n5Data/ezequiel23/a530ce83268de2a.n5?dataSetName=3340371272",
new String[]{"Dark", "Flour"},
new HashMap<Integer, String>(){{put(0, "ec"); put(1, "cytosol"); put(2, "Nucleus");}},
new double[][]{
{2.2250738585072014E-308, 2.2250738585072014E-308, 2.2250738585072014E-308, 3.52810806343277, 3.67238125879494,
0.37504412562962713, 0.1941393552238283},
{10.0, 10.0, 10.0, 10.0, 10.0, 10.0, 10.0}
},
new double[][]{
{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
{10.0, 10.0, 10.0, 6.471891936567229, 6.32761874120506, 9.624955874370375, 9.805860644776171}
},
new double[][]{
{0.0, 0.0, 0.0, 0.007280282509686344, 0.013997372023832945, 0.013997372023832947, 0.013997372023832949},
{9.999999999999494, 9.999999999999494, 9.999999999999494, 9.992719717489994, 9.986002627975711,
9.98600262797569, 9.986002627975722}
},
143300,
new double[]{0, 14891, 3697}
);

return new N5DataSetFile[]{frapSimulationResultsMasked, anns5DTIRFSimulation};
}

}
Original file line number Diff line number Diff line change
@@ -1,17 +1,21 @@
package org.vcell.N5;

import com.google.gson.internal.LinkedTreeMap;
import ij.ImagePlus;
import ij.io.Opener;
import ij.plugin.Duplicator;
import ij.plugin.ImageCalculator;
import ij.process.ImageProcessor;
import org.junit.*;
import org.junit.Assert;
import org.junit.Before;
import org.junit.Test;

import java.io.File;
import java.io.IOException;
import java.net.URISyntaxException;
import java.net.URL;
import java.util.ArrayList;
import java.util.Map;


/*
Expand Down Expand Up @@ -86,9 +90,13 @@ public void testS3AlphaInstanceLoadedIntoMemory() throws IOException {
SimResultsLoader simResultsLoader = new SimResultsLoader(n5DataSetFile.uri, "");
simResultsLoader.createS3Client();
ImagePlus imagePlus = simResultsLoader.getImgPlusFromN5File();
alphaStatsTest(new Duplicator().run(imagePlus), n5DataSetFile, stats.HISTMAX);
alphaStatsTest(new Duplicator().run(imagePlus), n5DataSetFile, stats.HISTMIN);
alphaStatsTest(new Duplicator().run(imagePlus), n5DataSetFile, stats.HISTAVERAGE);
ImagePlus inMemory = new Duplicator().run(imagePlus);
for (Object property : imagePlus.getProperties().keySet()){
inMemory.setProperty((String) property, imagePlus.getProperty((String) property));
}
alphaStatsTest(inMemory, n5DataSetFile, stats.HISTMAX);
alphaStatsTest(inMemory, n5DataSetFile, stats.HISTMIN);
alphaStatsTest(inMemory, n5DataSetFile, stats.HISTAVERAGE);
}
}

Expand All @@ -105,35 +113,100 @@ public void testUnits() throws IOException {
ImagePlus imagePlus = simResultsLoader.getImgPlusFromN5File();
double areaOfPixel = imagePlus.getCalibration().getX(1) * imagePlus.getCalibration().getY(1);
double totalArea = areaOfPixel * imagePlus.getWidth() * imagePlus.getHeight();

boolean threeD = imagePlus.getNSlices() > 1;
if (threeD){
areaOfPixel *= imagePlus.getCalibration().getZ(1);
totalArea *= imagePlus.getDimensions()[3];
}
Assert.assertEquals(n5DataSetFile.totalArea, totalArea, 0.0000001);

totalArea = 0;
imagePlus.setPosition(imagePlus.getNChannels(), 1, 1);
ImageProcessor imageProcessor = imagePlus.getProcessor();

PixelCalculations pixelCalculations = ((w, h, w1, h1) -> {
int inBounds = (0 <= w1) && (w1 < imagePlus.getWidth()) && (0 <= h1) && (h1 < imagePlus.getHeight()) ? 2 : 1;
if (inBounds == 2 && imageProcessor.getf(w1, h1) != imageProcessor.getf(w, h)){
boolean edgeOfDomain = (inBounds == 2) &&
(imageProcessor.getf(w1, h1) != imageProcessor.getf(w, h));
if (edgeOfDomain){
return inBounds;
}
return 1;
});
for (int h = 0; h < imagePlus.getHeight(); h++){
for (int w = 0; w < imagePlus.getWidth(); w++){
double sum = pixelCalculations.grabEdge(w, h, w, h -1) *

for(int domain: n5DataSetFile.mask.keySet()){
totalArea = 0;
for (int h = 0; h < imagePlus.getHeight(); h++){
for (int w = 0; w < imagePlus.getWidth(); w++){
double sum = pixelCalculations.grabEdge(w, h, w, h -1) *
pixelCalculations.grabEdge(w, h, w, h + 1) *
pixelCalculations.grabEdge(w, h,w + 1, h) *
pixelCalculations.grabEdge(w, h,w -1, h);
if (imageProcessor.getf(w, h) != 1 && sum > 1){
totalArea += (areaOfPixel / sum);
if (imageProcessor.getf(w,h) != domain && sum > 1 && domain != 0){
totalArea += (areaOfPixel / sum);
}
if (imageProcessor.getf(w, h) == domain && domain != 0){
totalArea += areaOfPixel;
}
}
}
Assert.assertEquals(n5DataSetFile.testDomainArea[domain], totalArea, n5DataSetFile.testDomainArea[domain] * 0.013);
}
}
}

private double alphaStatsThreeD(ImagePlus imagePlus, stats testType, int channel, int frame){
double experimentValue = Double.NaN;
for (int z = 1; z <= imagePlus.getNSlices(); z++){
imagePlus.setPosition(channel, z, frame);
double currentValue;
switch (testType){
case HISTMAX:
setImageMask(imagePlus);
currentValue = imagePlus.getStatistics().max;
if (currentValue == -Double.MAX_VALUE){
currentValue = 0;
}
if (imageProcessor.getf(w, h) == 1){
totalArea += areaOfPixel;
if (Double.isNaN(experimentValue) || experimentValue < currentValue){
experimentValue = currentValue;
}
break;
case HISTMIN:
setImageMask(imagePlus);
currentValue = imagePlus.getStatistics().min;
if (Double.isNaN(experimentValue) || experimentValue > currentValue){
experimentValue = currentValue;
}
break;
case HISTAVERAGE:
experimentValue = average3D(imagePlus);
}
if (testType.equals(stats.HISTAVERAGE)){
break;
}
}
return experimentValue;
}

private double average3D(ImagePlus imagePlus){
double total = 0;
double voxelVolume = imagePlus.getCalibration().pixelHeight * imagePlus.getCalibration().pixelWidth * imagePlus.getCalibration().pixelDepth;
double totalVolume = 0;
for (int z = 1; z < imagePlus.getNSlices(); z++){
imagePlus.setPosition(imagePlus.getChannel(), z, imagePlus.getFrame());
setImageMask(imagePlus);
ImageProcessor imageProcessor = imagePlus.getProcessor();
for (int x = 0; x < imagePlus.getWidth(); x++){
for (int y = 0; y < imagePlus.getHeight(); y++){
double pixelValue = imageProcessor.getValue(x, y);
if (!Double.isNaN(pixelValue)){
total += (pixelValue * voxelVolume);
totalVolume += voxelVolume;
}
}
}
Assert.assertEquals(n5DataSetFile.testDomainArea, totalArea, n5DataSetFile.testDomainArea * 0.013);

}
return (total / totalVolume );
}

public void alphaStatsTest(ImagePlus imagePlus, N5DataSetFile n5DataSetFile, stats testType){
Expand All @@ -148,25 +221,29 @@ public void alphaStatsTest(ImagePlus imagePlus, N5DataSetFile n5DataSetFile, sta
case HISTAVERAGE:
controlData = n5DataSetFile.histAverage;
}
for(int k = 0; k < controlData.length; k++){
for (int i = 0; i < controlData[k].length; i++){
imagePlus.setPosition(k + 1, 1, i + 1); //frame position seems to start at 1
boolean threeD = imagePlus.getNSlices() > 1;
for(int channel = 1; channel <= controlData.length; channel++){
for (int frame = 1; frame <= controlData[channel - 1].length; frame++){
imagePlus.setPosition(channel, 1, frame); //frame position seems to start at 1
double experimentalData = 0;
switch (testType){
case HISTMAX:
experimentalData = imagePlus.getStatistics().histMax;
experimentalData = threeD ? alphaStatsThreeD(imagePlus, testType, channel, frame) :
imagePlus.getStatistics().max;
break;
case HISTMIN:
setImageMask(imagePlus);
experimentalData = imagePlus.getStatistics().min;
experimentalData = threeD ? alphaStatsThreeD(imagePlus, testType, channel, frame) :
imagePlus.getStatistics().min;
break;
case HISTAVERAGE:
setImageMask(imagePlus);
experimentalData = imagePlus.getStatistics().mean;
experimentalData = threeD ? alphaStatsThreeD(imagePlus, testType, channel, frame) :
imagePlus.getStatistics().mean;
}
Assert.assertEquals("Channel: " + k + " Time: " + i + " Stat: " + testType +
"\n Experiment Value: " + experimentalData + " Control Value: " + controlData[k][i]
,0.0, experimentalData - controlData[k][i], 0.000001);
Assert.assertEquals("Channel: " + channel + " Time: " + frame + " Stat: " + testType +
"\n Experiment Value: " + experimentalData + " Control Value: " + controlData[channel - 1][frame - 1]
,0.0, experimentalData - controlData[channel - 1][frame - 1], 0.0001);
}
}
}
Expand Down

0 comments on commit b20a80e

Please sign in to comment.