Commit f654ff3d authored by Lucy McNeill's avatar Lucy McNeill
Browse files

count_foci has individual mask routines. More modular and closer to...

count_foci has individual mask routines. More modular and closer to bioconductor standards. Mask routines can be called with different settings for exceptional cases e.g. after removing XY, low background and low signal
parent b3156eb9
No preview for this file type
......@@ -47,6 +47,7 @@
#' @param artificial_amp_factor Amplification of foci channel, for annotation only.
#' @param strand_amp multiplication of strand channel to make masks
#' @param min_foci minimum pixel area for a foci. Depends on your dpi etc. Defaults to 4
#' @param disc_size size of disc for local background calculation in dna channel
#' @examples demo_path = paste0(system.file("extdata",package = "synapsis"))
#' foci_counts <- count_foci(demo_path,offset_factor = 3, brush_size = 3,
#' brush_sigma = 3, annotation = "on",stage = "pachytene")
......@@ -54,11 +55,10 @@
#' @return foci count per cell
count_foci <- function(img_path, stage = "none", offset_px = 0.2, offset_factor = 2, brush_size = 3, brush_sigma = 3, foci_norm = 0.01, annotation = "off",channel2_string = "SYCP3", channel1_string = "MLH3",file_ext = "jpeg", KO_str = "--",WT_str = "++",KO_out = "-/-", WT_out = "+/+", watershed_stop = "off", watershed_radius = 1, watershed_tol = 0.05, crowded_foci = TRUE, artificial_amp_factor = 1, strand_amp = 2, min_foci =-1)
count_foci <- function(img_path, stage = "none", offset_px = 0.2, offset_factor = 2, brush_size = 3, brush_sigma = 3, foci_norm = 0.01, annotation = "off",channel2_string = "SYCP3", channel1_string = "MLH3",file_ext = "jpeg", KO_str = "--",WT_str = "++",KO_out = "-/-", WT_out = "+/+", watershed_stop = "off", watershed_radius = 1, watershed_tol = 0.05, crowded_foci = TRUE, artificial_amp_factor = 1, strand_amp = 2, min_foci =-1, disc_size = 51)
{
cell_count <- 0
image_count <-0
foci_counts <- 0
antibody1_store <- 0
antibody2_store <- 0
if(stage == "pachytene"){
......@@ -98,74 +98,32 @@ count_foci <- function(img_path, stage = "none", offset_px = 0.2, offset_factor
antibody1_store <- 0
antibody2_store <- 0
cell_count <- cell_count + 1
new_img<-img_orig
disc <- makeBrush(21, "disc")
disc <- disc / sum(disc)
localBackground <- filter2(new_img, disc)
offset <- offset_px
if(stage == "pachytene"){
thresh_crop <- (new_img - localBackground > offset)
}
else{
thresh_crop <- new_img > offset
}
strands <- bwlabel(thresh_crop)
#### call mask strand channel
strands <- make_strand_mask(offset_px, stage, img_orig, disc_size)
color_img_strands<- colorLabels(strands, normalize = TRUE)
num_strands <- computeFeatures.shape(strands)
num_strands <- data.frame(num_strands)
foci_mask_crop <- img_orig_foci
### call mask foci foci channel
bg <- mean(img_orig_foci)
orig_mean <- mean(img_orig_foci)
mean_factor <- foci_norm/orig_mean
img_orig_foci <- img_orig_foci*mean_factor
#### normalise the foci image
offset <- offset_factor*bg
if(crowded_foci == TRUE){
foci_th <- foci_mask_crop > bg + offset
}
else{
### smooth it
img_tmp_contrast <- foci_mask_crop
w <- makeBrush(size = brush_size, shape = 'gaussian', sigma = brush_sigma)
img_flo <- filter2(img_tmp_contrast, w)
## smooth foci channel
foci_th <- img_flo > bg + offset
}
foci_label <- bwlabel(foci_th)
foci_label <- make_foci_mask(offset_factor,bg,crowded_foci,img_orig_foci,brush_size,brush_sigma)
foci_label <- channel(foci_label, "grey")
num_strands <- computeFeatures.shape(strands)
num_strands <- data.frame(num_strands)
##### print properties of the images
if(watershed_stop != "off"){
coincident_foci <- bwlabel(foci_label*strands)
}
else{
coincident_foci <- watershed(bwlabel(foci_th*strands)*as.matrix(img_orig_foci),tolerance=watershed_tol, ext=watershed_radius)
}
coincident_df <- data.frame(computeFeatures.shape(coincident_foci))
if(annotation == "on"){
print(coincident_df)
}
coincident_df <- coincident_df[coincident_df$s.area > min_foci,]
### multiply strands by foci_label
if(annotation == "on"){
annotate_foci_counting(img_file,cell_count,new_img,img_orig_foci,artificial_amp_factor,foci_mask_crop,strands,coincident_foci, foci_label)
}
overlap_no <- table(coincident_foci)
foci_per_cell <- length(overlap_no)
foci_per_cell <- nrow(coincident_df)
foci_candidates <- computeFeatures.shape(foci_label)
foci_candidates <- data.frame(foci_candidates)
foci_areas <- foci_candidates$s.area
coincident_foci <- get_overlap_mask(strands, foci_label, watershed_stop, img_orig_foci, watershed_radius, watershed_tol)
foci_per_cell <- get_foci_per_cell(img_file,offset_px,stage,strands,watershed_stop,foci_label,annotation, min_foci,cell_count,img_orig, img_orig_foci, artificial_amp_factor,coincident_foci)
#### end foci counting function
if(annotation=="on"){
cat("\n which counts this many foci:",foci_per_cell, sep = " ")
}
image_mat <- as.matrix(foci_mask_crop)
image_mat <- as.matrix(img_orig_foci)
image_mat <- image_mat[image_mat > 1e-06]
mean_ratio <- median(image_mat)/mean(image_mat)
skew <- (median(image_mat)-mean(image_mat))/sd(image_mat)
### look at properties of the foci.
foci_candidates <- computeFeatures.shape(foci_label)
foci_candidates <- data.frame(foci_candidates)
foci_areas <- foci_candidates$s.area
### look at properties of the overlap foci.
##
overlap_candidates <- computeFeatures.shape(coincident_foci)
overlap_candidates <- data.frame(overlap_candidates)
overlap <- overlap_candidates$s.area
......@@ -180,14 +138,25 @@ count_foci <- function(img_path, stage = "none", offset_px = 0.2, offset_factor
}
}
percent_px <- sum(overlap)/sum(foci_areas)
### if percent_px too small:
if(annotation == "on"){
cat("\n percentage of foci channel coincident:", percent_px*100, sep = " ")
# count_low_contrast_image()
tryCatch({
if(percent_px < 0.3){
cat("\n file name", file_foci, "coincides", percent_px*100,
"percent. So decided to reduce the local background thresholding
requirements.", sep = " ")
foci_candidates <- make_foci_mask(0.5*offset_factor,0.5*bg,crowded_foci,img_orig_foci,brush_size,brush_sigma)
### if it's zero percent... might be genuinely zero. Ignore?
### if it's non zero, largish: probably noisy. No distinction between
### high background/ high signal vs low bg low signal..
}
},
error = function(e) {
#what should be done in case of exception?
cat("\n percent on was a NaN. Suspect no/low signal")
}
)
df_cells <- append_data_frame(WT_str,KO_str,WT_out,KO_out,img_file,foci_areas,df_cells,cell_count,stage,foci_per_cell,image_mat,percent_px,alone_foci)
}
}
colnames(df_cells) <- df_cols
......@@ -202,28 +171,27 @@ count_foci <- function(img_path, stage = "none", offset_px = 0.2, offset_factor
#'
#' @param img_file cell's file name
#' @param cell_count unique cell counter
#' @param new_img cropped strand channel
#' @param img_orig original strand crop
#' @param img_orig_foci cropped foci channel
#' @param artificial_amp_factor amplification factor
#' @param foci_mask_crop black white mask of foci channel
#' @param strands black white mask of strand channel
#' @param coincident_foci mask of overlap between strand and foci channel
#' @param foci_label black and white mask of foci channel
#'
#'
annotate_foci_counting <- function(img_file,cell_count,new_img,img_orig_foci,artificial_amp_factor,foci_mask_crop,strands,coincident_foci, foci_label){
annotate_foci_counting <- function(img_file,cell_count,img_orig,img_orig_foci,artificial_amp_factor,strands,coincident_foci, foci_label){
cat("\n at file",img_file, sep = " ")
cat("\n cell counter is", cell_count, sep= " ")
cat("\n original images")
plot(new_img)
plot(img_orig)
plot(img_orig_foci)
ch1 <-channel(new_img,"grey")
ch2 <- channel(artificial_amp_factor*foci_mask_crop,"grey")
ch1 <-channel(img_orig,"grey")
ch2 <- channel(artificial_amp_factor*img_orig_foci,"grey")
bluered <- rgbImage(ch1, ch2, 0*ch1)
cat("\n displaying resulting foci count plots. Overlay two channels:")
plot(rgbImage(ch1,ch2,0*new_img))
plot(rgbImage(ch1,ch2,0*img_orig))
cat("\n displaying resulting masks. Overlay two masks:")
plot(rgbImage(strands,foci_label,0*new_img))
plot(rgbImage(strands,foci_label,0*img_orig))
cat("\n coincident foci:")
plot(colorLabels(coincident_foci))
cat("\n two channels, only coincident foci")
......@@ -281,4 +249,149 @@ append_data_frame <- function(WT_str,KO_str,WT_out,KO_out,img_file,foci_areas,df
}
### annotate when the flag has been raised that
### annotate when the flag has been raised that the percent of foci on is too small.
## Nan case indicates zero, so should move to a separate case of dim foci possibly.
## Change the thresholding to longer be local background (but still smoothed)
## perhaps base the new threshold value on the mean pixel value... rather then difference compared to background?
### then make a test about whether there was an improvement....
#### annonate when a flag has been raised about the XY stuff:
### if sd(foci) is >10* typical value, then you have an XY blob there. Remove it from the mask?
#' remove_XY
#'
#' applies new row to data frame
#'
#' @param mask current foci mask
#'
#'
remove_XY <- function(){
### if the standard deviation was too big
### loop over and delete the giant blobs.
### after this function is run, coincident foci needs to be computed again.
}
#' make_foci_mask
#'
#' creates foci mask for foci channel crop
#'
#' @param offset_factor Pixel value offset used in thresholding of foci channel
#' @param bg background value- currently just mean pixel value of whole image
#' @param crowded_foci TRUE or FALSE, defaults to FALSE. Set to TRUE if you
#' have foci > 100 or so.
#' @param img_orig_foci cropped foci channel
#' @param brush_size size of brush to smooth the foci channel. Should be small
#' to avoid erasing foci.
#' @param brush_sigma sigma for Gaussian smooth of foci channel. Should be
#' small to avoid erasing foci.
#'
#'
#'
make_foci_mask <- function(offset_factor,bg,crowded_foci,img_orig_foci,brush_size,brush_sigma){
foci_mask_crop <- img_orig_foci
offset <- offset_factor*bg
if(crowded_foci == TRUE){
foci_th <- foci_mask_crop > bg + offset
}
else{
### smooth it
img_tmp_contrast <- foci_mask_crop
w <- makeBrush(size = brush_size, shape = 'gaussian', sigma = brush_sigma)
img_flo <- filter2(img_tmp_contrast, w)
## smooth foci channel
foci_th <- img_flo > bg + offset
}
foci_label <- bwlabel(foci_th)
return(foci_label)
# here check whether there is a huge blob. Remove it
# remove_XY()
#
}
#' make_strand_mask
#'
#' creates strand mask for strand channel crop
#'
#' @param offset_px, Pixel value offset used in thresholding of dna channel
#' @param stage meitoic stage, currently pachytene or not.
#' @param img_orig original strand crop
#' @param disc_size size of disc for local background calculation in dna channel
#'
#'
make_strand_mask <- function(offset_px, stage, img_orig, disc_size){
new_img<-img_orig
disc <- makeBrush(disc_size, "disc")
disc <- disc / sum(disc)
localBackground <- filter2(new_img, disc)
offset <- offset_px
if(stage == "pachytene"){
thresh_crop <- (new_img - localBackground > offset)
}
else{
thresh_crop <- new_img > offset
}
strands <- bwlabel(thresh_crop)
return(strands)
}
#' get_overlap_mask
#'
#' creates mask for coincident foci
#'
#' @param strands black white mask of strand channel
#' @param watershed_stop Turn off default watershed method with "off"
#' @param foci_label black and white mask of foci channel
#' @param img_orig_foci cropped foci channel
#' @param watershed_radius Radius (ext variable) in watershed method used
#' in foci channel. Defaults to 1 (small)
#' @param watershed_tol Intensity tolerance for watershed method. Defaults to 0.05.
#'
get_overlap_mask<- function(strands, foci_label, watershed_stop, img_orig_foci, watershed_radius, watershed_tol){
num_strands <- computeFeatures.shape(strands)
num_strands <- data.frame(num_strands)
if(watershed_stop != "off"){
coincident_foci <- bwlabel(foci_label*strands)
}
else{
foci_th <- channel(foci_label,"grey")
coincident_foci <- watershed(bwlabel(foci_th*strands)*as.matrix(img_orig_foci),tolerance=watershed_tol, ext=watershed_radius)
}
return(coincident_foci)
}
#' get_foci_per_cell
#'
#' creates mask for coincident foci
#'
#' @param img_file cell's file name
#' @param offset_px, Pixel value offset used in thresholding of dna channel
#' @param stage meitoic stage, currently pachytene or not.
#' @param strands black white mask of strand channel
#' @param watershed_stop Turn off default watershed method with "off"
#' @param foci_label black and white mask of foci channel
#' @param annotation, Choice to output pipeline choices (recommended to knit)
#' @param min_foci minimum pixel area for a foci. Depends on your dpi etc. Defaults to 4
#' @param cell_count unique cell counter
#' @param img_orig original strand crop
#' @param img_orig_foci cropped foci channel
#' @param artificial_amp_factor amplification factor
#' @param coincident_foci mask of coincident foci
#'
#'
#'
get_foci_per_cell <- function(img_file,offset_px,stage,strands,watershed_stop,foci_label, annotation, min_foci, cell_count, img_orig, img_orig_foci, artificial_amp_factor, coincident_foci){
coincident_df <- data.frame(computeFeatures.shape(coincident_foci))
if(annotation == "on"){
print(coincident_df)
}
coincident_df <- coincident_df[coincident_df$s.area > min_foci,]
### multiply strands by foci_label
if(annotation == "on"){
annotate_foci_counting(img_file,cell_count,img_orig,img_orig_foci,artificial_amp_factor,strands,coincident_foci, foci_label)
}
foci_per_cell <- nrow(coincident_df)
return(foci_per_cell)
}
......@@ -7,10 +7,9 @@
annotate_foci_counting(
img_file,
cell_count,
new_img,
img_orig,
img_orig_foci,
artificial_amp_factor,
foci_mask_crop,
strands,
coincident_foci,
foci_label
......@@ -21,14 +20,12 @@ annotate_foci_counting(
\item{cell_count}{unique cell counter}
\item{new_img}{cropped strand channel}
\item{img_orig}{original strand crop}
\item{img_orig_foci}{cropped foci channel}
\item{artificial_amp_factor}{amplification factor}
\item{foci_mask_crop}{black white mask of foci channel}
\item{strands}{black white mask of strand channel}
\item{coincident_foci}{mask of overlap between strand and foci channel}
......
......@@ -26,7 +26,8 @@ count_foci(
crowded_foci = TRUE,
artificial_amp_factor = 1,
strand_amp = 2,
min_foci = -1
min_foci = -1,
disc_size = 51
)
}
\arguments{
......@@ -85,6 +86,8 @@ have foci > 100 or so.}
\item{strand_amp}{multiplication of strand channel to make masks}
\item{min_foci}{minimum pixel area for a foci. Depends on your dpi etc. Defaults to 4}
\item{disc_size}{size of disc for local background calculation in dna channel}
}
\value{
foci count per cell
......
% Generated by roxygen2: do not edit by hand
% Please edit documentation in R/count_foci.R
\name{get_foci_per_cell}
\alias{get_foci_per_cell}
\title{get_foci_per_cell}
\usage{
get_foci_per_cell(
img_file,
offset_px,
stage,
strands,
watershed_stop,
foci_label,
annotation,
min_foci,
cell_count,
img_orig,
img_orig_foci,
artificial_amp_factor,
coincident_foci
)
}
\arguments{
\item{img_file}{cell's file name}
\item{offset_px, }{Pixel value offset used in thresholding of dna channel}
\item{stage}{meitoic stage, currently pachytene or not.}
\item{strands}{black white mask of strand channel}
\item{watershed_stop}{Turn off default watershed method with "off"}
\item{foci_label}{black and white mask of foci channel}
\item{annotation, }{Choice to output pipeline choices (recommended to knit)}
\item{min_foci}{minimum pixel area for a foci. Depends on your dpi etc. Defaults to 4}
\item{cell_count}{unique cell counter}
\item{img_orig}{original strand crop}
\item{img_orig_foci}{cropped foci channel}
\item{artificial_amp_factor}{amplification factor}
\item{coincident_foci}{mask of coincident foci}
}
\description{
creates mask for coincident foci
}
% Generated by roxygen2: do not edit by hand
% Please edit documentation in R/count_foci.R
\name{get_overlap_mask}
\alias{get_overlap_mask}
\title{get_overlap_mask}
\usage{
get_overlap_mask(
strands,
foci_label,
watershed_stop,
img_orig_foci,
watershed_radius,
watershed_tol
)
}
\arguments{
\item{strands}{black white mask of strand channel}
\item{foci_label}{black and white mask of foci channel}
\item{watershed_stop}{Turn off default watershed method with "off"}
\item{img_orig_foci}{cropped foci channel}
\item{watershed_radius}{Radius (ext variable) in watershed method used
in foci channel. Defaults to 1 (small)}
\item{watershed_tol}{Intensity tolerance for watershed method. Defaults to 0.05.}
}
\description{
creates mask for coincident foci
}
% Generated by roxygen2: do not edit by hand
% Please edit documentation in R/count_foci.R
\name{make_foci_mask}
\alias{make_foci_mask}
\title{make_foci_mask}
\usage{
make_foci_mask(
offset_factor,
bg,
crowded_foci,
img_orig_foci,
brush_size,
brush_sigma
)
}
\arguments{
\item{offset_factor}{Pixel value offset used in thresholding of foci channel}
\item{bg}{background value- currently just mean pixel value of whole image}
\item{crowded_foci}{TRUE or FALSE, defaults to FALSE. Set to TRUE if you
have foci > 100 or so.}
\item{img_orig_foci}{cropped foci channel}
\item{brush_size}{size of brush to smooth the foci channel. Should be small
to avoid erasing foci.}
\item{brush_sigma}{sigma for Gaussian smooth of foci channel. Should be
small to avoid erasing foci.}
}
\description{
creates foci mask for foci channel crop
}
% Generated by roxygen2: do not edit by hand
% Please edit documentation in R/count_foci.R
\name{make_strand_mask}
\alias{make_strand_mask}
\title{make_strand_mask}
\usage{
make_strand_mask(offset_px, stage, img_orig, disc_size)
}
\arguments{
\item{offset_px, }{Pixel value offset used in thresholding of dna channel}
\item{stage}{meitoic stage, currently pachytene or not.}
\item{img_orig}{original strand crop}
\item{disc_size}{size of disc for local background calculation in dna channel}
}
\description{
creates strand mask for strand channel crop
}
% Generated by roxygen2: do not edit by hand
% Please edit documentation in R/count_foci.R
\name{remove_XY}
\alias{remove_XY}
\title{remove_XY}
\usage{
remove_XY()
}
\arguments{
\item{mask}{current foci mask}
}
\description{
applies new row to data frame
}
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment