In TKoscik/nifti.qc: Quality control metrics and reports for NIFTI data
rm(list=ls())
gc()
knitr::opts_chunk$set(echo=FALSE, warning=FALSE, message=FALSE, error=FALSE)
# Debug inputs ----
researcher <- "/Shared/nopoulos"
project <- "sca_pilot"
subject <- "sub-231"
session <- "ses-328zk16wb6"
site <- "site-00201"
space <- "HCPMNI2009c"
template <- "800um"
img.source <- "acq-sagMPRAGEPROMO_T1w"
scratch.dir <- "/Shared/koscikt_scratch/nifti_qc_scratch"
nimg.root <- "/Shared/nopoulos/nimg_core"
# ----
img.mod <- unlist(strsplit(img.source, "_"))
img.mod <- img.mod[length(img.mod)]
library(nifti.qc)
library(nifti.io)
library(nifti.draw)
library(R.utils)
library(ggplot2)
library(viridis)
library(reshape2)
library(gridExtra)
library(grid)
library(raster)
library(moments)
library(OpenImageR)
library(kableExtra)
Anatomical Preprocessing Report
INC, r format(Sys.time(), "%Y-%m-%d %H:%M:%S-%z")
Researcher: /Shared/research_lab
Project: project_name
Participant: sub-231
Session: ses-328zk16wb6
Raw Image: /Shared/researchlab/project_name/nifti/sub-231/ses-328zk16wb6/anat/sub-231_ses-328zk16wb6_site-00201_acq-sagMPRAGEPROMO_T1w.nii.gz
Processed Native: /Shared/researchlab/project_name/derivatives/anat/native/sub-231_ses-328zk16wb6_site-00201_T1w.nii.gz
Processed Normalized: /Shared/researchlab/project_name/derivatives/anat/reg_HCPMNI2009c_800um/sub-231_ses-328zk16wb6_site-00201_reg-HCPMNI2009c+800um_T1w.nii.gz
# gather needed images
source.file <- paste0(researcher, "/", project, "/nifti/", subject, "/", session, "/anat/", subject, "_", session, "_", site, "_", img.source, ".nii.gz")
img.rigid <- paste0(scratch.dir, "/", subject, "_", session, "_rigid.nii")
img.native <- paste0(scratch.dir, "/", subject, "_", session, "_native.nii")
img.template <- paste0(scratch.dir, "/template.nii")
mask.air <- paste0(scratch.dir, "/", subject, "_", session, "_mask-air.nii")
mask.brain <- paste0(scratch.dir, "/", subject, "_", session, "_mask-brain.nii")
mask.brain.reg <- paste0(scratch.dir, "/", subject, "_", session, "_mask-brain-reg.nii")
seg.label <- paste0(scratch.dir, "/", subject, "_", session, "_seg-label.nii")
prob.csf <- paste0(scratch.dir, "/", subject, "_", session, "_seg-CSF.nii")
prob.gm <- paste0(scratch.dir, "/", subject, "_", session, "_seg-GM.nii")
prob.wm <- paste0(scratch.dir, "/", subject, "_", session, "_seg-WM.nii")
img.noise <- paste0(scratch.dir, "/", subject, "_", session, "_noise.nii")
img.biast1t2 <- paste0(scratch.dir, "/", subject, "_", session, "_biasFieldT1T2.nii")
img.biasn4 <- paste0(scratch.dir, "/", subject, "_", session, "_biasFieldN4.nii")
img.art.rigid <- paste0(scratch.dir, "/", subject, "_", session, "_artRigid.nii")
img.art.native <- paste0(scratch.dir, "/", subject, "_", session, "_artNative.nii")
gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/prep/", subject, "/", session, "/", subject, "_", session, "_", site, "_", img.source, "_prep-rigid.nii.gz"), destname=img.rigid, remove=FALSE, overwrite=TRUE)
gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/native/", subject, "_", session, "_", site, "_", img.mod, ".nii.gz"), destname=img.native, remove=FALSE, overwrite=TRUE)
gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/mask/", subject, "_", session, "_", site, "_mask-air.nii.gz"), destname=mask.air, remove=FALSE, overwrite=TRUE)
gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/mask/", subject, "_", session, "_", site, "_mask-brain.nii.gz"), destname=mask.brain, remove=FALSE, overwrite=TRUE)
gunzip(filename=paste0(nimg.root, "/templates/", space, "/", template, "/", space, "_", template, "_", img.mod, ".nii.gz"), destname=img.template, remove=FALSE, overwrite=TRUE)
gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/reg_", space, "_", template, "/", subject, "_", session, "_", site, "_reg-", space, "+", template, "_T1w_brain.nii.gz"), destname=mask.brain.reg, remove=FALSE, overwrite=TRUE)
gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/segmentation/", subject, "_", session, "_", site, "_seg-label.nii.gz"), destname=seg.label, remove=FALSE, overwrite=TRUE)
gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/segmentation/", subject, "_", session, "_", site, "_seg-CSF.nii.gz"), destname=prob.csf, remove=FALSE, overwrite=TRUE)
gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/segmentation/", subject, "_", session, "_", site, "_seg-GM.nii.gz"), destname=prob.gm, remove=FALSE, overwrite=TRUE)
gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/segmentation/", subject, "_", session, "_", site, "_seg-WM.nii.gz"), destname=prob.wm, remove=FALSE, overwrite=TRUE)
gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/prep/", subject, "/", session, "/", subject, "_", session, "_", site, "_", img.source, "_prep-noise.nii.gz"), destname=img.noise, remove=FALSE, overwrite=TRUE)
gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/prep/", subject, "/", session, "/", subject, "_", session, "_", site, "_prep-biasFieldT1T2.nii.gz"), destname=img.biast1t2, remove=FALSE, overwrite=TRUE)
gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/prep/", subject, "/", session, "/", subject, "_", session, "_", site, "_", img.mod, "_prep-biasFieldN4.nii.gz"), destname=img.biasn4, remove=FALSE, overwrite=TRUE)
wbf <- data.frame(
Metric = c("SNR", "SNR.D", "CNR", "QI.1", "QI.2", "CJV",
"FWHM.average", "FWHM.x", "FWHM.y", "FWHM.z",
"EFC", "FBER", "WM.to.MAX"),
Raw.Value=numeric(13),
Clean.Value=numeric(13))
var.num <- 0
var.num <- var.num + 1
wbf$Raw.Value[var.num] <- nii.qc.snr(
img.nii=img.rigid, img.vol=1,
mask.nii=mask.brain, mask.vol=1, mask.dir="gt", mask.thresh=0)
wbf$Clean.Value[var.num] <- nii.qc.snr(
img.nii=img.native, img.vol=1,
mask.nii=mask.brain, mask.vol=1, mask.dir="gt", mask.thresh=0)
var.num <- var.num + 1
wbf$Raw.Value[var.num] <- nii.qc.snrd(
img.nii=img.rigid, img.vol=1,
mask.nii=mask.brain, mask.vol=1, mask.dir="gt", mask.thresh=0,
air.nii=mask.air, air.vol=1, air.dir="eq", air.thresh=1)
wbf$Clean.Value[var.num] <- nii.qc.snrd(
img.nii=img.native, img.vol=1,
mask.nii=mask.brain, mask.vol=1, mask.dir="gt", mask.thresh=0,
air.nii=mask.air, air.vol=1, air.dir="eq", air.thresh=1)
var.num <- var.num + 1
wbf$Raw.Value[var.num] <- nii.qc.cnr(img.nii=img.rigid, img.vol=1,
gm.nii=seg.label, gm.vol=1, gm.dir="eq", gm.thresh=2,
wm.nii=seg.label, wm.vol=1, wm.dir="eq", wm.thresh=3,
air.nii=mask.air, air.vol=1, air.dir="eq", air.thresh=1)
wbf$Clean.Value[var.num] <- nii.qc.cnr(img.nii=img.native, img.vol=1,
gm.nii=seg.label, gm.vol=1, gm.dir="eq", gm.thresh=2,
wm.nii=seg.label, wm.vol=1, wm.dir="eq", wm.thresh=3,
air.nii=mask.air, air.vol=1, air.dir="eq", air.thresh=1)
invisible(nii.qc.artefacts(
img.nii=img.rigid, img.vol=1L,
air.nii=mask.air, air.vol=1L, air.dir="eq", air.thresh=1,
save.dir=scratch.dir, file.name=paste0(subject, "_", session, "_artRigid.nii")))
qi <- nii.qc.mortamet(
img.nii=img.rigid, img.vol=1L,
air.nii=mask.air, air.vol=1L, air.dir="eq", air.thresh=1,
art.nii=img.art.rigid, art.vol=1L, art.dir="eq", art.thresh=1)
var.num <- var.num + 1
wbf$Raw.Value[var.num] <- qi$qi1
wbf$Clean.Value[var.num] <- NA
var.num <- var.num + 1
wbf$Raw.Value[var.num] <- qi$qi2
wbf$Clean.Value[var.num] <- NA
var.num <- var.num + 1
wbf$Raw.Value[var.num] <- nii.qc.cjv(
img.nii=img.rigid, img.vol=1L,
gm.nii=seg.label, gm.vol=1L, gm.dir="eq", gm.thresh=2,
wm.nii=seg.label, wm.vol=1L, wm.dir="eq", wm.thresh=3)
wbf$Clean.Value[var.num] <- nii.qc.cjv(
img.nii=img.native, img.vol=1L,
gm.nii=seg.label, gm.vol=1L, gm.dir="eq", gm.thresh=2,
wm.nii=seg.label, wm.vol=1L, wm.dir="eq", wm.thresh=3)
fwhm.rigid <- nii.qc.fwhm(img.nii=img.rigid, mask.nii=mask.brain, save.dir=scratch.dir, file.name="rigid.fwhm.txt")
fwhm.native <- nii.qc.fwhm(img.nii=img.native, mask.nii=mask.brain, save.dir=scratch.dir, file.name="native.fwhm.txt")
var.num <- var.num + 1
wbf$Raw.Value[var.num] <- fwhm.rigid$average
wbf$Clean.Value[var.num] <- fwhm.native$average
var.num <- var.num + 1
wbf$Raw.Value[var.num] <- fwhm.rigid$x
wbf$Clean.Value[var.num] <- fwhm.native$x
var.num <- var.num + 1
wbf$Raw.Value[var.num] <- fwhm.rigid$y
wbf$Clean.Value[var.num] <- fwhm.native$y
var.num <- var.num + 1
wbf$Raw.Value[var.num] <- fwhm.rigid$z
wbf$Clean.Value[var.num] <- fwhm.native$z
var.num <- var.num + 1
wbf$Raw.Value[var.num] <- nii.qc.efc(img.nii=img.rigid, img.vol=1L)
wbf$Clean.Value[var.num] <- nii.qc.efc(img.nii=img.native, img.vol=1L)
var.num <- var.num + 1
wbf$Raw.Value[var.num] <- nii.qc.fber(img.nii=img.rigid, img.vol=1L, brain.mask=mask.brain, brain.vol=1L, brain.dir="gt", brain.thresh = 0)
wbf$Clean.Value[var.num] <- nii.qc.fber(img.nii=img.native, img.vol=1L, brain.mask=mask.brain, brain.vol=1L, brain.dir="gt", brain.thresh = 0)
var.num <- var.num + 1
wbf$Raw.Value[var.num] <- nii.qc.wm2max (img.nii=img.rigid, img.vol=1L, wm.nii=seg.label, wm.vol=1L, wm.dir="eq", wm.thresh=3)
wbf$Clean.Value[var.num] <- nii.qc.wm2max (img.nii=img.native, img.vol=1L, wm.nii=seg.label, wm.vol=1L, wm.dir="eq", wm.thresh=3)
segf <- data.frame(
Metric = c("Ratio.to.Whole.Brain", "rPVe", "SNR", "SNR.D",
"Mean", "SD", "Median", "MAD", "Skewness", "Kurtosis",
"quantile.05", "quantile.95"),
CSF=numeric(12),
GM=numeric(12),
WM=numeric(12))
segf$CSF[1] <- nii.qc.volfrac(
whole.nii=seg.label, whole.vol=1L, whole.dir="gt", whole.thresh=1,
part.nii=seg.label, part.vol=1L, part.dir="eq", part.thresh=1)
segf$GM[1] <- nii.qc.volfrac(
whole.nii=seg.label, whole.vol=1L, whole.dir="gt", whole.thresh=1,
part.nii=seg.label, part.vol=1L, part.dir="eq", part.thresh=2)
segf$WM[1] <- nii.qc.volfrac(
whole.nii=seg.label, whole.vol=1L, whole.dir="gt", whole.thresh=1,
part.nii=seg.label, part.vol=1L, part.dir="eq", part.thresh=3)
segf$CSF[2] <- nii.qc.rpve(tissue.prob=prob.csf, tissue.vol=1L)
segf$GM[2] <- nii.qc.rpve(tissue.prob=prob.gm, tissue.vol=1L)
segf$WM[2] <- nii.qc.rpve(tissue.prob=prob.wm, tissue.vol=1L)
segf$CSF[3] <- nii.qc.snr(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=1)
segf$GM[3] <- nii.qc.snr(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=2)
segf$WM[3] <- nii.qc.snr(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=3)
segf$CSF[4] <- nii.qc.snrd(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=1, air.nii=mask.air, air.vol=1L, air.dir="eq", air.thresh=1)
segf$GM[4] <- nii.qc.snrd(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=2, air.nii=mask.air, air.vol=1L, air.dir="eq", air.thresh=1)
segf$WM[4] <- nii.qc.snrd(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=3, air.nii=mask.air, air.vol=1L, air.dir="eq", air.thresh=1)
segf$CSF[5:12] <- as.numeric(nii.qc.descriptive(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=1)$stats)
segf$GM[5:12] <- as.numeric(nii.qc.descriptive(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=2)$stats)
segf$WM[5:12] <- as.numeric(nii.qc.descriptive(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=3)$stats)
inc.norms.wb <- read.csv(paste0(nimg.root, "/qc_norms/inc_anat_wb.csv"), as.is=TRUE)
if (source.file %in% inc.norms.wb$input.image) {
inc.norms.wb[which(source.file == inc.norms.wb$input.image)[1], 3:15] <- as.matrix(t(wbf)[2,])
inc.norms.wb[which(source.file == inc.norms.wb$input.image)[2], 3:15] <- as.matrix(t(wbf)[3,])
} else {
tf <- data.frame(input.image = rep(source.file, 2),
Value = c("Raw", "Clean"),
SNR=as.numeric(wbf[1,2:3]),
SNR.D=as.numeric(wbf[2,2:3]),
CNR=as.numeric(wbf[3,2:3]),
QI.1=as.numeric(wbf[4,2:3]),
QI.2=as.numeric(wbf[5,2:3]),
CJV=as.numeric(wbf[6,2:3]),
FWHM.average=as.numeric(wbf[7,2:3]),
FWHM.x=as.numeric(wbf[8,2:3]),
FWHM.y=as.numeric(wbf[9,2:3]),
FWHM.z=as.numeric(wbf[10,2:3]),
EFC=as.numeric(wbf[11,2:3]),
FBER=as.numeric(wbf[12,2:3]),
WM.to.MAX=as.numeric(wbf[13,2:3]))
inc.norms.wb <- rbind(inc.norms.wb, tf)
}
write.table(x=inc.norms.wb, file=paste0(nimg.root, "/qc_norms/inc_anat_wb.csv"),
row.names = FALSE, col.names=TRUE, quote = FALSE, sep = ",")
norm.vals <- read.csv(paste0(nimg.root, "/qc_norms/inc_anat_wb.csv"), as.is=TRUE)
plotf <- wbf
colnames(plotf) <- c("Metric", "Raw", "Clean")
plotf <- melt(plotf, id.vars="Metric")
plotf$lower <- numeric(nrow(plotf))*NA
plotf$upper <- numeric(nrow(plotf))*NA
plotf$norm <- numeric(nrow(plotf))
for (i in 1:nrow(plotf)) {
if (!is.na(plotf$value[i])) {
norm.m <- mean(norm.vals[(norm.vals$Value==plotf$variable[i]), which(colnames(norm.vals)==plotf$Metric[i])])
norm.sd <- sd(norm.vals[(norm.vals$Value==plotf$variable[i]), which(colnames(norm.vals)==plotf$Metric[i])])
iqr <- IQR(norm.vals[(norm.vals$Value==plotf$variable[i]), which(colnames(norm.vals)==plotf$Metric[i])])
plotf$value[i] <- (plotf$value[i] - norm.m) /norm.sd
plotf$lower[i] <- (-iqr) / norm.sd
plotf$upper[i] <- (iqr) / norm.sd
}
}
plotf$ok <- factor(((plotf$value > plotf$lower) & (plotf$value < plotf$upper))*1, levels=c(0,1,NA), labels=c("Check", "OK"))
plotf$Metric <- factor(plotf$Metric, levels=c("WM.to.MAX", "FBER", "EFC", "FWHM.z", "FWHM.y", "FWHM.x", "FWHM.average", "CJV", "QI.2", "QI.1", "CNR", "SNR.D", "SNR"))
if (!("Check" %in% plotf$ok)) {
ok.colors <- "#64a8ff"
} else if (!("OK" %in% plotf$ok)) {
ok.colors <- "#ff0000"
} else {
ok.colors <- c("#ff0000", "#64a8ff")
}
plot1 <- ggplot(plotf, aes(x=Metric)) +
theme_bw() +
geom_pointrange(aes(y=norm, ymin=lower, ymax=upper, group=variable, color=variable),
position=position_dodge(width=0.3), shape=124, size=1,
show.legend = FALSE) +
scale_shape_manual(values=c(21,23)) +
scale_color_manual(values=c("#a8a8a8", "#000000"), guide="none") +
scale_fill_manual(values=ok.colors, guide="none") +
geom_point(aes(y=value, shape=variable, group=variable, fill=ok),
size=2, position=position_dodge(width=0.3)) +
coord_flip() +
labs(subtitle="Whole Brain", x="", y="Z-Score") +
theme(legend.title = element_blank(),
legend.position = "bottom",
legend.direction = "horizontal",
legend.background = element_blank(),
legend.text = element_text(size=8),
plot.title = element_blank(),
plot.subtitle = element_text(size=10),
axis.title = element_text(size=10),
axis.text.x = element_text(size=8),
axis.text.y = element_text(size=8),
plot.margin=margin(0,0,0,0, "null"),
legend.margin = margin(0,0,0,0, "null"))
inc.norms.seg <- read.csv(paste0(nimg.root, "/qc_norms/inc_anat_seg.csv"), as.is=TRUE)
if (source.file %in% inc.norms.seg$input.image) {
inc.norms.seg[which(source.file == inc.norms.seg$input.image)[1], 3:14] <- as.matrix(t(segf)[2,])
inc.norms.seg[which(source.file == inc.norms.seg$input.image)[2], 3:14] <- as.matrix(t(segf)[3,])
inc.norms.seg[which(source.file == inc.norms.seg$input.image)[3], 3:14] <- as.matrix(t(segf)[4,])
} else {
tf <- data.frame(input.image = rep(source.file, 3),
Tissue = c("CSF", "GM", "WM"),
Ratio.to.Whole.Brain = as.numeric(segf[1, 2:4]),
rPVe = as.numeric(segf[2, 2:4]),
SNR = as.numeric(segf[3, 2:4]),
SNR.d = as.numeric(segf[4, 2:4]),
Mean = as.numeric(segf[5, 2:4]),
SD = as.numeric(segf[6, 2:4]),
Median = as.numeric(segf[7, 2:4]),
MAD = as.numeric(segf[8, 2:4]),
Skewness = as.numeric(segf[9, 2:4]),
Kurtosis = as.numeric(segf[10, 2:4]),
quantile.05 = as.numeric(segf[11, 2:4]),
quantile.95 = as.numeric(segf[12, 2:4]))
inc.norms.seg <- rbind(inc.norms.seg, tf)
}
write.table(x=inc.norms.seg, file=paste0(nimg.root, "/qc_norms/inc_anat_seg.csv"),
row.names = FALSE, col.names=TRUE, quote = FALSE, sep = ",")
norm.vals <- read.csv(paste0(nimg.root, "/qc_norms/inc_anat_seg.csv"), as.is=TRUE)
plotf <- segf
plotf <- melt(plotf, id.vars="Metric")
plotf$lower <- numeric(nrow(plotf))*NA
plotf$upper <- numeric(nrow(plotf))*NA
plotf$norm <- numeric(nrow(plotf))
for (i in 1:nrow(plotf)) {
if (!is.na(plotf$value[i])) {
norm.m <- mean(norm.vals[(norm.vals$Tissue==plotf$variable[i]), which(colnames(norm.vals)==plotf$Metric[i])])
norm.sd <- sd(norm.vals[(norm.vals$Tissue==plotf$variable[i]), which(colnames(norm.vals)==plotf$Metric[i])])
iqr <- IQR(norm.vals[(norm.vals$Tissue==plotf$variable[i]), which(colnames(norm.vals)==plotf$Metric[i])])
plotf$value[i] <- (plotf$value[i] - norm.m) /norm.sd
plotf$lower[i] <- (-iqr) / norm.sd
plotf$upper[i] <- (iqr) / norm.sd
}
}
plotf$ok <- factor(((plotf$value > plotf$lower) & (plotf$value < plotf$upper))*1, levels=c(0,1, NA), labels=c("Check", "OK"))
plotf$Metric <- factor(plotf$Metric, levels=c("quantile.95", "quantile.05", "Kurtosis", "Skewness", "MAD", "Median", "SD", "Mean", "rPVe", "Ratio.to.Whole.Brain", "SNR.D", "SNR"))
if (!("Check" %in% plotf$ok)) {
ok.colors <- "#64a8ff"
} else if (!("OK" %in% plotf$ok)) {
ok.colors <- "#ff0000"
} else {
ok.colors <- c("#ff0000", "#64a8ff")
}
plot2 <- ggplot(plotf, aes(x=Metric)) +
theme_bw() +
geom_pointrange(aes(y=norm, ymin=lower, ymax=upper, group=variable, color=variable),
position=position_dodge(width=0.3), shape=124, size=1,
show.legend = FALSE) +
scale_shape_manual(values=c(21,23,22)) +
scale_color_manual(values=c("#000000", "#646464", "#a8a8a8"), guide="none") +
scale_fill_manual(values=ok.colors, guide="none") +
coord_flip() +
geom_point(aes(y=value, shape=variable, group=variable, fill=ok),
size=2, position=position_dodge(width=0.3)) +
labs(subtitle="by Tissue Class", x="", y="Z-Score") +
theme(legend.title = element_blank(),
legend.position = "bottom",
legend.direction = "horizontal",
legend.background = element_blank(),
legend.text = element_text(size=8),
plot.title = element_blank(),
plot.subtitle = element_text(size=10),
axis.title = element_text(size=10),
axis.text.x = element_text(size=8),
axis.text.y = element_text(size=8),
plot.margin=margin(0,0,0,0, "null"),
legend.margin = margin(0,0,0,0, "null"))
grid.arrange(textGrob("Quality Metrics", rot = 90), plot1, plot2, ncol=, widths=c(0.1,1,1))
volume <- read.nii.volume(img.rigid, 1L)
coords <- c(round(dim(volume)[1]/5*2), round(dim(volume)[2:3]/2))
slice1 <- volume[coords[1], , ]
slice2 <- volume[ , coords[2], ]
slice3 <- volume[ , , coords[3]]
pixdims <- unlist(nii.hdr(img.rigid, "pixdim"))[2:4]
slice1 <- resizeImage(slice1, width=dim(slice1)[1]*pixdims[2], height=dim(slice1)[2]*pixdims[3])
slice2 <- resizeImage(slice2, width=dim(slice2)[1]*pixdims[1], height=dim(slice2)[2]*pixdims[3])
slice3 <- resizeImage(slice3, width=dim(slice3)[1]*pixdims[1], height=dim(slice3)[2]*pixdims[2])
slices <- rbind(slice1, slice2, slice3)
slices <- melt(slices)
# plot1 <-
ggplot(slices, aes(x=Var1, y=Var2, fill=value)) +
theme_bw() +
coord_equal(expand=FALSE, clip="off") +
scale_fill_gradient(low="#000000", high="#ffffff") +
geom_raster() +
annotate("label", x=-Inf, y=Inf, label="Unprocessed",
hjust=0, vjust=1, color="#000000", fill="#ffffff",
label.r=unit(0, "null")) +
theme(plot.title = element_blank(),
legend.position="none",
legend.title = element_blank(),
legend.text = element_text(size=8, margin=margin(1,0,0,0,"null")),
axis.title=element_blank(),
axis.text=element_blank(),
axis.ticks=element_blank(),
plot.subtitle = element_text(size=10, margin = margin(0,0,0,0,"null")),
plot.background=element_blank(),
panel.background=element_rect(color="#000000"),
panel.grid=element_blank(),
panel.border=element_blank(),
panel.spacing.x=unit(c(0,0,0,0),"null"),
panel.spacing.y=unit(c(0,0,0,0),"null"),
plot.margin=margin(0,0,0,0, "null"),
legend.margin=margin(0,0,0,0, "null"),
panel.spacing = margin(0,0,0,0, "null"))
rm(list=c("volume", "coords", "slice1", "slice2", "slice3", "slices", "pixdims"))
volume <- read.nii.volume(img.native, 1L)
coords <- c(round(dim(volume)[1]/5*2), round(dim(volume)[2:3]/2))
slice1 <- volume[coords[1], , ]
slice2 <- volume[ , coords[2], ]
slice3 <- volume[ , , coords[3]]
pixdims <- unlist(nii.hdr(img.rigid, "pixdim"))[2:4]
slice1 <- resizeImage(slice1, width=dim(slice1)[1]*pixdims[2], height=dim(slice1)[2]*pixdims[3])
slice2 <- resizeImage(slice2, width=dim(slice2)[1]*pixdims[1], height=dim(slice2)[2]*pixdims[3])
slice3 <- resizeImage(slice3, width=dim(slice3)[1]*pixdims[1], height=dim(slice3)[2]*pixdims[2])
slices <- rbind(slice1, slice2, slice3)
slices <- melt(slices)
# plot2 <-
ggplot(slices, aes(x=Var1, y=Var2, fill=value)) +
theme_bw() +
coord_equal(expand=FALSE, clip="off") +
scale_fill_gradient(low="#000000", high="#ffffff") +
geom_raster() +
annotate("label", x=-Inf, y=Inf, label="Cleaned",
hjust=0, vjust=1, color="#000000", fill="#ffffff",
label.r=unit(0, "null")) +
theme(plot.title = element_blank(),
legend.position="none",
legend.title = element_blank(),
legend.text = element_blank(),
axis.title=element_blank(),
axis.text=element_blank(),
axis.ticks=element_blank(),
plot.background=element_blank(),
panel.background=element_blank(),
panel.grid=element_blank(),
panel.border=element_blank(),
panel.spacing.x=unit(c(0,0,0,0),"null"),
panel.spacing.y=unit(c(0,0,0,0),"null"),
plot.margin=margin(0,0,0,0, "null"),
legend.margin=margin(0,0,0,0, "null"),
panel.spacing = margin(0,0,0,0, "null"))
rm(list=c("volume", "coords", "slice1", "slice2", "slice3", "slices", "pixdims"))
wzxhzdk:11
wzxhzdk:12
wzxhzdk:13
wzxhzdk:14
wzxhzdk:15
wzxhzdk:16
wzxhzdk:17
wzxhzdk:18
wzxhzdk:19
wzxhzdk:20
TKoscik/nifti.qc documentation built on May 6, 2019, 4:30 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
rm(list=ls()) gc() knitr::opts_chunk$set(echo=FALSE, warning=FALSE, message=FALSE, error=FALSE) # Debug inputs ---- researcher <- "/Shared/nopoulos" project <- "sca_pilot" subject <- "sub-231" session <- "ses-328zk16wb6" site <- "site-00201" space <- "HCPMNI2009c" template <- "800um" img.source <- "acq-sagMPRAGEPROMO_T1w" scratch.dir <- "/Shared/koscikt_scratch/nifti_qc_scratch" nimg.root <- "/Shared/nopoulos/nimg_core" # ---- img.mod <- unlist(strsplit(img.source, "_")) img.mod <- img.mod[length(img.mod)]
library(nifti.qc) library(nifti.io) library(nifti.draw) library(R.utils) library(ggplot2) library(viridis) library(reshape2) library(gridExtra) library(grid) library(raster) library(moments) library(OpenImageR) library(kableExtra)
Anatomical Preprocessing Report
INC, r format(Sys.time(), "%Y-%m-%d %H:%M:%S-%z")
Researcher: /Shared/research_lab
Project: project_name
Participant: sub-231
Session: ses-328zk16wb6
Raw Image: /Shared/researchlab/project_name/nifti/sub-231/ses-328zk16wb6/anat/sub-231_ses-328zk16wb6_site-00201_acq-sagMPRAGEPROMO_T1w.nii.gz
Processed Native: /Shared/researchlab/project_name/derivatives/anat/native/sub-231_ses-328zk16wb6_site-00201_T1w.nii.gz
Processed Normalized: /Shared/researchlab/project_name/derivatives/anat/reg_HCPMNI2009c_800um/sub-231_ses-328zk16wb6_site-00201_reg-HCPMNI2009c+800um_T1w.nii.gz
# gather needed images source.file <- paste0(researcher, "/", project, "/nifti/", subject, "/", session, "/anat/", subject, "_", session, "_", site, "_", img.source, ".nii.gz") img.rigid <- paste0(scratch.dir, "/", subject, "_", session, "_rigid.nii") img.native <- paste0(scratch.dir, "/", subject, "_", session, "_native.nii") img.template <- paste0(scratch.dir, "/template.nii") mask.air <- paste0(scratch.dir, "/", subject, "_", session, "_mask-air.nii") mask.brain <- paste0(scratch.dir, "/", subject, "_", session, "_mask-brain.nii") mask.brain.reg <- paste0(scratch.dir, "/", subject, "_", session, "_mask-brain-reg.nii") seg.label <- paste0(scratch.dir, "/", subject, "_", session, "_seg-label.nii") prob.csf <- paste0(scratch.dir, "/", subject, "_", session, "_seg-CSF.nii") prob.gm <- paste0(scratch.dir, "/", subject, "_", session, "_seg-GM.nii") prob.wm <- paste0(scratch.dir, "/", subject, "_", session, "_seg-WM.nii") img.noise <- paste0(scratch.dir, "/", subject, "_", session, "_noise.nii") img.biast1t2 <- paste0(scratch.dir, "/", subject, "_", session, "_biasFieldT1T2.nii") img.biasn4 <- paste0(scratch.dir, "/", subject, "_", session, "_biasFieldN4.nii") img.art.rigid <- paste0(scratch.dir, "/", subject, "_", session, "_artRigid.nii") img.art.native <- paste0(scratch.dir, "/", subject, "_", session, "_artNative.nii")
gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/prep/", subject, "/", session, "/", subject, "_", session, "_", site, "_", img.source, "_prep-rigid.nii.gz"), destname=img.rigid, remove=FALSE, overwrite=TRUE) gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/native/", subject, "_", session, "_", site, "_", img.mod, ".nii.gz"), destname=img.native, remove=FALSE, overwrite=TRUE) gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/mask/", subject, "_", session, "_", site, "_mask-air.nii.gz"), destname=mask.air, remove=FALSE, overwrite=TRUE) gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/mask/", subject, "_", session, "_", site, "_mask-brain.nii.gz"), destname=mask.brain, remove=FALSE, overwrite=TRUE) gunzip(filename=paste0(nimg.root, "/templates/", space, "/", template, "/", space, "_", template, "_", img.mod, ".nii.gz"), destname=img.template, remove=FALSE, overwrite=TRUE) gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/reg_", space, "_", template, "/", subject, "_", session, "_", site, "_reg-", space, "+", template, "_T1w_brain.nii.gz"), destname=mask.brain.reg, remove=FALSE, overwrite=TRUE) gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/segmentation/", subject, "_", session, "_", site, "_seg-label.nii.gz"), destname=seg.label, remove=FALSE, overwrite=TRUE) gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/segmentation/", subject, "_", session, "_", site, "_seg-CSF.nii.gz"), destname=prob.csf, remove=FALSE, overwrite=TRUE) gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/segmentation/", subject, "_", session, "_", site, "_seg-GM.nii.gz"), destname=prob.gm, remove=FALSE, overwrite=TRUE) gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/segmentation/", subject, "_", session, "_", site, "_seg-WM.nii.gz"), destname=prob.wm, remove=FALSE, overwrite=TRUE) gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/prep/", subject, "/", session, "/", subject, "_", session, "_", site, "_", img.source, "_prep-noise.nii.gz"), destname=img.noise, remove=FALSE, overwrite=TRUE) gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/prep/", subject, "/", session, "/", subject, "_", session, "_", site, "_prep-biasFieldT1T2.nii.gz"), destname=img.biast1t2, remove=FALSE, overwrite=TRUE) gunzip(filename=paste0(researcher, "/", project, "/derivatives/anat/prep/", subject, "/", session, "/", subject, "_", session, "_", site, "_", img.mod, "_prep-biasFieldN4.nii.gz"), destname=img.biasn4, remove=FALSE, overwrite=TRUE)
wbf <- data.frame( Metric = c("SNR", "SNR.D", "CNR", "QI.1", "QI.2", "CJV", "FWHM.average", "FWHM.x", "FWHM.y", "FWHM.z", "EFC", "FBER", "WM.to.MAX"), Raw.Value=numeric(13), Clean.Value=numeric(13)) var.num <- 0 var.num <- var.num + 1 wbf$Raw.Value[var.num] <- nii.qc.snr( img.nii=img.rigid, img.vol=1, mask.nii=mask.brain, mask.vol=1, mask.dir="gt", mask.thresh=0) wbf$Clean.Value[var.num] <- nii.qc.snr( img.nii=img.native, img.vol=1, mask.nii=mask.brain, mask.vol=1, mask.dir="gt", mask.thresh=0) var.num <- var.num + 1 wbf$Raw.Value[var.num] <- nii.qc.snrd( img.nii=img.rigid, img.vol=1, mask.nii=mask.brain, mask.vol=1, mask.dir="gt", mask.thresh=0, air.nii=mask.air, air.vol=1, air.dir="eq", air.thresh=1) wbf$Clean.Value[var.num] <- nii.qc.snrd( img.nii=img.native, img.vol=1, mask.nii=mask.brain, mask.vol=1, mask.dir="gt", mask.thresh=0, air.nii=mask.air, air.vol=1, air.dir="eq", air.thresh=1) var.num <- var.num + 1 wbf$Raw.Value[var.num] <- nii.qc.cnr(img.nii=img.rigid, img.vol=1, gm.nii=seg.label, gm.vol=1, gm.dir="eq", gm.thresh=2, wm.nii=seg.label, wm.vol=1, wm.dir="eq", wm.thresh=3, air.nii=mask.air, air.vol=1, air.dir="eq", air.thresh=1) wbf$Clean.Value[var.num] <- nii.qc.cnr(img.nii=img.native, img.vol=1, gm.nii=seg.label, gm.vol=1, gm.dir="eq", gm.thresh=2, wm.nii=seg.label, wm.vol=1, wm.dir="eq", wm.thresh=3, air.nii=mask.air, air.vol=1, air.dir="eq", air.thresh=1) invisible(nii.qc.artefacts( img.nii=img.rigid, img.vol=1L, air.nii=mask.air, air.vol=1L, air.dir="eq", air.thresh=1, save.dir=scratch.dir, file.name=paste0(subject, "_", session, "_artRigid.nii"))) qi <- nii.qc.mortamet( img.nii=img.rigid, img.vol=1L, air.nii=mask.air, air.vol=1L, air.dir="eq", air.thresh=1, art.nii=img.art.rigid, art.vol=1L, art.dir="eq", art.thresh=1) var.num <- var.num + 1 wbf$Raw.Value[var.num] <- qi$qi1 wbf$Clean.Value[var.num] <- NA var.num <- var.num + 1 wbf$Raw.Value[var.num] <- qi$qi2 wbf$Clean.Value[var.num] <- NA var.num <- var.num + 1 wbf$Raw.Value[var.num] <- nii.qc.cjv( img.nii=img.rigid, img.vol=1L, gm.nii=seg.label, gm.vol=1L, gm.dir="eq", gm.thresh=2, wm.nii=seg.label, wm.vol=1L, wm.dir="eq", wm.thresh=3) wbf$Clean.Value[var.num] <- nii.qc.cjv( img.nii=img.native, img.vol=1L, gm.nii=seg.label, gm.vol=1L, gm.dir="eq", gm.thresh=2, wm.nii=seg.label, wm.vol=1L, wm.dir="eq", wm.thresh=3) fwhm.rigid <- nii.qc.fwhm(img.nii=img.rigid, mask.nii=mask.brain, save.dir=scratch.dir, file.name="rigid.fwhm.txt") fwhm.native <- nii.qc.fwhm(img.nii=img.native, mask.nii=mask.brain, save.dir=scratch.dir, file.name="native.fwhm.txt") var.num <- var.num + 1 wbf$Raw.Value[var.num] <- fwhm.rigid$average wbf$Clean.Value[var.num] <- fwhm.native$average var.num <- var.num + 1 wbf$Raw.Value[var.num] <- fwhm.rigid$x wbf$Clean.Value[var.num] <- fwhm.native$x var.num <- var.num + 1 wbf$Raw.Value[var.num] <- fwhm.rigid$y wbf$Clean.Value[var.num] <- fwhm.native$y var.num <- var.num + 1 wbf$Raw.Value[var.num] <- fwhm.rigid$z wbf$Clean.Value[var.num] <- fwhm.native$z var.num <- var.num + 1 wbf$Raw.Value[var.num] <- nii.qc.efc(img.nii=img.rigid, img.vol=1L) wbf$Clean.Value[var.num] <- nii.qc.efc(img.nii=img.native, img.vol=1L) var.num <- var.num + 1 wbf$Raw.Value[var.num] <- nii.qc.fber(img.nii=img.rigid, img.vol=1L, brain.mask=mask.brain, brain.vol=1L, brain.dir="gt", brain.thresh = 0) wbf$Clean.Value[var.num] <- nii.qc.fber(img.nii=img.native, img.vol=1L, brain.mask=mask.brain, brain.vol=1L, brain.dir="gt", brain.thresh = 0) var.num <- var.num + 1 wbf$Raw.Value[var.num] <- nii.qc.wm2max (img.nii=img.rigid, img.vol=1L, wm.nii=seg.label, wm.vol=1L, wm.dir="eq", wm.thresh=3) wbf$Clean.Value[var.num] <- nii.qc.wm2max (img.nii=img.native, img.vol=1L, wm.nii=seg.label, wm.vol=1L, wm.dir="eq", wm.thresh=3)
segf <- data.frame( Metric = c("Ratio.to.Whole.Brain", "rPVe", "SNR", "SNR.D", "Mean", "SD", "Median", "MAD", "Skewness", "Kurtosis", "quantile.05", "quantile.95"), CSF=numeric(12), GM=numeric(12), WM=numeric(12)) segf$CSF[1] <- nii.qc.volfrac( whole.nii=seg.label, whole.vol=1L, whole.dir="gt", whole.thresh=1, part.nii=seg.label, part.vol=1L, part.dir="eq", part.thresh=1) segf$GM[1] <- nii.qc.volfrac( whole.nii=seg.label, whole.vol=1L, whole.dir="gt", whole.thresh=1, part.nii=seg.label, part.vol=1L, part.dir="eq", part.thresh=2) segf$WM[1] <- nii.qc.volfrac( whole.nii=seg.label, whole.vol=1L, whole.dir="gt", whole.thresh=1, part.nii=seg.label, part.vol=1L, part.dir="eq", part.thresh=3) segf$CSF[2] <- nii.qc.rpve(tissue.prob=prob.csf, tissue.vol=1L) segf$GM[2] <- nii.qc.rpve(tissue.prob=prob.gm, tissue.vol=1L) segf$WM[2] <- nii.qc.rpve(tissue.prob=prob.wm, tissue.vol=1L) segf$CSF[3] <- nii.qc.snr(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=1) segf$GM[3] <- nii.qc.snr(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=2) segf$WM[3] <- nii.qc.snr(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=3) segf$CSF[4] <- nii.qc.snrd(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=1, air.nii=mask.air, air.vol=1L, air.dir="eq", air.thresh=1) segf$GM[4] <- nii.qc.snrd(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=2, air.nii=mask.air, air.vol=1L, air.dir="eq", air.thresh=1) segf$WM[4] <- nii.qc.snrd(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=3, air.nii=mask.air, air.vol=1L, air.dir="eq", air.thresh=1) segf$CSF[5:12] <- as.numeric(nii.qc.descriptive(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=1)$stats) segf$GM[5:12] <- as.numeric(nii.qc.descriptive(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=2)$stats) segf$WM[5:12] <- as.numeric(nii.qc.descriptive(img.nii=img.native, img.vol=1L, mask.nii=seg.label, mask.vol=1L, mask.dir="eq", mask.thresh=3)$stats)
inc.norms.wb <- read.csv(paste0(nimg.root, "/qc_norms/inc_anat_wb.csv"), as.is=TRUE) if (source.file %in% inc.norms.wb$input.image) { inc.norms.wb[which(source.file == inc.norms.wb$input.image)[1], 3:15] <- as.matrix(t(wbf)[2,]) inc.norms.wb[which(source.file == inc.norms.wb$input.image)[2], 3:15] <- as.matrix(t(wbf)[3,]) } else { tf <- data.frame(input.image = rep(source.file, 2), Value = c("Raw", "Clean"), SNR=as.numeric(wbf[1,2:3]), SNR.D=as.numeric(wbf[2,2:3]), CNR=as.numeric(wbf[3,2:3]), QI.1=as.numeric(wbf[4,2:3]), QI.2=as.numeric(wbf[5,2:3]), CJV=as.numeric(wbf[6,2:3]), FWHM.average=as.numeric(wbf[7,2:3]), FWHM.x=as.numeric(wbf[8,2:3]), FWHM.y=as.numeric(wbf[9,2:3]), FWHM.z=as.numeric(wbf[10,2:3]), EFC=as.numeric(wbf[11,2:3]), FBER=as.numeric(wbf[12,2:3]), WM.to.MAX=as.numeric(wbf[13,2:3])) inc.norms.wb <- rbind(inc.norms.wb, tf) } write.table(x=inc.norms.wb, file=paste0(nimg.root, "/qc_norms/inc_anat_wb.csv"), row.names = FALSE, col.names=TRUE, quote = FALSE, sep = ",") norm.vals <- read.csv(paste0(nimg.root, "/qc_norms/inc_anat_wb.csv"), as.is=TRUE) plotf <- wbf colnames(plotf) <- c("Metric", "Raw", "Clean") plotf <- melt(plotf, id.vars="Metric") plotf$lower <- numeric(nrow(plotf))*NA plotf$upper <- numeric(nrow(plotf))*NA plotf$norm <- numeric(nrow(plotf)) for (i in 1:nrow(plotf)) { if (!is.na(plotf$value[i])) { norm.m <- mean(norm.vals[(norm.vals$Value==plotf$variable[i]), which(colnames(norm.vals)==plotf$Metric[i])]) norm.sd <- sd(norm.vals[(norm.vals$Value==plotf$variable[i]), which(colnames(norm.vals)==plotf$Metric[i])]) iqr <- IQR(norm.vals[(norm.vals$Value==plotf$variable[i]), which(colnames(norm.vals)==plotf$Metric[i])]) plotf$value[i] <- (plotf$value[i] - norm.m) /norm.sd plotf$lower[i] <- (-iqr) / norm.sd plotf$upper[i] <- (iqr) / norm.sd } } plotf$ok <- factor(((plotf$value > plotf$lower) & (plotf$value < plotf$upper))*1, levels=c(0,1,NA), labels=c("Check", "OK")) plotf$Metric <- factor(plotf$Metric, levels=c("WM.to.MAX", "FBER", "EFC", "FWHM.z", "FWHM.y", "FWHM.x", "FWHM.average", "CJV", "QI.2", "QI.1", "CNR", "SNR.D", "SNR")) if (!("Check" %in% plotf$ok)) { ok.colors <- "#64a8ff" } else if (!("OK" %in% plotf$ok)) { ok.colors <- "#ff0000" } else { ok.colors <- c("#ff0000", "#64a8ff") } plot1 <- ggplot(plotf, aes(x=Metric)) + theme_bw() + geom_pointrange(aes(y=norm, ymin=lower, ymax=upper, group=variable, color=variable), position=position_dodge(width=0.3), shape=124, size=1, show.legend = FALSE) + scale_shape_manual(values=c(21,23)) + scale_color_manual(values=c("#a8a8a8", "#000000"), guide="none") + scale_fill_manual(values=ok.colors, guide="none") + geom_point(aes(y=value, shape=variable, group=variable, fill=ok), size=2, position=position_dodge(width=0.3)) + coord_flip() + labs(subtitle="Whole Brain", x="", y="Z-Score") + theme(legend.title = element_blank(), legend.position = "bottom", legend.direction = "horizontal", legend.background = element_blank(), legend.text = element_text(size=8), plot.title = element_blank(), plot.subtitle = element_text(size=10), axis.title = element_text(size=10), axis.text.x = element_text(size=8), axis.text.y = element_text(size=8), plot.margin=margin(0,0,0,0, "null"), legend.margin = margin(0,0,0,0, "null"))
inc.norms.seg <- read.csv(paste0(nimg.root, "/qc_norms/inc_anat_seg.csv"), as.is=TRUE) if (source.file %in% inc.norms.seg$input.image) { inc.norms.seg[which(source.file == inc.norms.seg$input.image)[1], 3:14] <- as.matrix(t(segf)[2,]) inc.norms.seg[which(source.file == inc.norms.seg$input.image)[2], 3:14] <- as.matrix(t(segf)[3,]) inc.norms.seg[which(source.file == inc.norms.seg$input.image)[3], 3:14] <- as.matrix(t(segf)[4,]) } else { tf <- data.frame(input.image = rep(source.file, 3), Tissue = c("CSF", "GM", "WM"), Ratio.to.Whole.Brain = as.numeric(segf[1, 2:4]), rPVe = as.numeric(segf[2, 2:4]), SNR = as.numeric(segf[3, 2:4]), SNR.d = as.numeric(segf[4, 2:4]), Mean = as.numeric(segf[5, 2:4]), SD = as.numeric(segf[6, 2:4]), Median = as.numeric(segf[7, 2:4]), MAD = as.numeric(segf[8, 2:4]), Skewness = as.numeric(segf[9, 2:4]), Kurtosis = as.numeric(segf[10, 2:4]), quantile.05 = as.numeric(segf[11, 2:4]), quantile.95 = as.numeric(segf[12, 2:4])) inc.norms.seg <- rbind(inc.norms.seg, tf) } write.table(x=inc.norms.seg, file=paste0(nimg.root, "/qc_norms/inc_anat_seg.csv"), row.names = FALSE, col.names=TRUE, quote = FALSE, sep = ",") norm.vals <- read.csv(paste0(nimg.root, "/qc_norms/inc_anat_seg.csv"), as.is=TRUE) plotf <- segf plotf <- melt(plotf, id.vars="Metric") plotf$lower <- numeric(nrow(plotf))*NA plotf$upper <- numeric(nrow(plotf))*NA plotf$norm <- numeric(nrow(plotf)) for (i in 1:nrow(plotf)) { if (!is.na(plotf$value[i])) { norm.m <- mean(norm.vals[(norm.vals$Tissue==plotf$variable[i]), which(colnames(norm.vals)==plotf$Metric[i])]) norm.sd <- sd(norm.vals[(norm.vals$Tissue==plotf$variable[i]), which(colnames(norm.vals)==plotf$Metric[i])]) iqr <- IQR(norm.vals[(norm.vals$Tissue==plotf$variable[i]), which(colnames(norm.vals)==plotf$Metric[i])]) plotf$value[i] <- (plotf$value[i] - norm.m) /norm.sd plotf$lower[i] <- (-iqr) / norm.sd plotf$upper[i] <- (iqr) / norm.sd } } plotf$ok <- factor(((plotf$value > plotf$lower) & (plotf$value < plotf$upper))*1, levels=c(0,1, NA), labels=c("Check", "OK")) plotf$Metric <- factor(plotf$Metric, levels=c("quantile.95", "quantile.05", "Kurtosis", "Skewness", "MAD", "Median", "SD", "Mean", "rPVe", "Ratio.to.Whole.Brain", "SNR.D", "SNR")) if (!("Check" %in% plotf$ok)) { ok.colors <- "#64a8ff" } else if (!("OK" %in% plotf$ok)) { ok.colors <- "#ff0000" } else { ok.colors <- c("#ff0000", "#64a8ff") } plot2 <- ggplot(plotf, aes(x=Metric)) + theme_bw() + geom_pointrange(aes(y=norm, ymin=lower, ymax=upper, group=variable, color=variable), position=position_dodge(width=0.3), shape=124, size=1, show.legend = FALSE) + scale_shape_manual(values=c(21,23,22)) + scale_color_manual(values=c("#000000", "#646464", "#a8a8a8"), guide="none") + scale_fill_manual(values=ok.colors, guide="none") + coord_flip() + geom_point(aes(y=value, shape=variable, group=variable, fill=ok), size=2, position=position_dodge(width=0.3)) + labs(subtitle="by Tissue Class", x="", y="Z-Score") + theme(legend.title = element_blank(), legend.position = "bottom", legend.direction = "horizontal", legend.background = element_blank(), legend.text = element_text(size=8), plot.title = element_blank(), plot.subtitle = element_text(size=10), axis.title = element_text(size=10), axis.text.x = element_text(size=8), axis.text.y = element_text(size=8), plot.margin=margin(0,0,0,0, "null"), legend.margin = margin(0,0,0,0, "null"))
grid.arrange(textGrob("Quality Metrics", rot = 90), plot1, plot2, ncol=, widths=c(0.1,1,1))
volume <- read.nii.volume(img.rigid, 1L) coords <- c(round(dim(volume)[1]/5*2), round(dim(volume)[2:3]/2)) slice1 <- volume[coords[1], , ] slice2 <- volume[ , coords[2], ] slice3 <- volume[ , , coords[3]] pixdims <- unlist(nii.hdr(img.rigid, "pixdim"))[2:4] slice1 <- resizeImage(slice1, width=dim(slice1)[1]*pixdims[2], height=dim(slice1)[2]*pixdims[3]) slice2 <- resizeImage(slice2, width=dim(slice2)[1]*pixdims[1], height=dim(slice2)[2]*pixdims[3]) slice3 <- resizeImage(slice3, width=dim(slice3)[1]*pixdims[1], height=dim(slice3)[2]*pixdims[2]) slices <- rbind(slice1, slice2, slice3) slices <- melt(slices) # plot1 <- ggplot(slices, aes(x=Var1, y=Var2, fill=value)) + theme_bw() + coord_equal(expand=FALSE, clip="off") + scale_fill_gradient(low="#000000", high="#ffffff") + geom_raster() + annotate("label", x=-Inf, y=Inf, label="Unprocessed", hjust=0, vjust=1, color="#000000", fill="#ffffff", label.r=unit(0, "null")) + theme(plot.title = element_blank(), legend.position="none", legend.title = element_blank(), legend.text = element_text(size=8, margin=margin(1,0,0,0,"null")), axis.title=element_blank(), axis.text=element_blank(), axis.ticks=element_blank(), plot.subtitle = element_text(size=10, margin = margin(0,0,0,0,"null")), plot.background=element_blank(), panel.background=element_rect(color="#000000"), panel.grid=element_blank(), panel.border=element_blank(), panel.spacing.x=unit(c(0,0,0,0),"null"), panel.spacing.y=unit(c(0,0,0,0),"null"), plot.margin=margin(0,0,0,0, "null"), legend.margin=margin(0,0,0,0, "null"), panel.spacing = margin(0,0,0,0, "null")) rm(list=c("volume", "coords", "slice1", "slice2", "slice3", "slices", "pixdims"))
volume <- read.nii.volume(img.native, 1L) coords <- c(round(dim(volume)[1]/5*2), round(dim(volume)[2:3]/2)) slice1 <- volume[coords[1], , ] slice2 <- volume[ , coords[2], ] slice3 <- volume[ , , coords[3]] pixdims <- unlist(nii.hdr(img.rigid, "pixdim"))[2:4] slice1 <- resizeImage(slice1, width=dim(slice1)[1]*pixdims[2], height=dim(slice1)[2]*pixdims[3]) slice2 <- resizeImage(slice2, width=dim(slice2)[1]*pixdims[1], height=dim(slice2)[2]*pixdims[3]) slice3 <- resizeImage(slice3, width=dim(slice3)[1]*pixdims[1], height=dim(slice3)[2]*pixdims[2]) slices <- rbind(slice1, slice2, slice3) slices <- melt(slices) # plot2 <- ggplot(slices, aes(x=Var1, y=Var2, fill=value)) + theme_bw() + coord_equal(expand=FALSE, clip="off") + scale_fill_gradient(low="#000000", high="#ffffff") + geom_raster() + annotate("label", x=-Inf, y=Inf, label="Cleaned", hjust=0, vjust=1, color="#000000", fill="#ffffff", label.r=unit(0, "null")) + theme(plot.title = element_blank(), legend.position="none", legend.title = element_blank(), legend.text = element_blank(), axis.title=element_blank(), axis.text=element_blank(), axis.ticks=element_blank(), plot.background=element_blank(), panel.background=element_blank(), panel.grid=element_blank(), panel.border=element_blank(), panel.spacing.x=unit(c(0,0,0,0),"null"), panel.spacing.y=unit(c(0,0,0,0),"null"), plot.margin=margin(0,0,0,0, "null"), legend.margin=margin(0,0,0,0, "null"), panel.spacing = margin(0,0,0,0, "null")) rm(list=c("volume", "coords", "slice1", "slice2", "slice3", "slices", "pixdims"))
wzxhzdk:11
wzxhzdk:12
wzxhzdk:13
wzxhzdk:14 wzxhzdk:15 wzxhzdk:16 wzxhzdk:17 wzxhzdk:18 wzxhzdk:19
wzxhzdk:20
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.