qra/python/extract_suv_from_mask.py

# -*- coding: utf-8 -*-
"""
Created on Thu Sep 29 14:51:09 2022

@author: katja
"""

import numpy
import os
import pandas as pd
import nibabel
import csv
import matplotlib.pyplot as plt

def organ_percentile(img, mask, level, p):
#ORGAN_PERCENTILE - computes suv_x the pth percentile of the distribution of
#   image intensity values in img from within some ROI mask
#
#   Inputs:
#       img - image. ndarray
#       mask - ROI mask. Must be same dimension as img. Must be binary ndarray
#       p - percentile. Between 0-100
#
#   Outputs:
#       suv_x - percentile of distribution. Defined by:
#
#                 suv_x
#           p/100 = ∫ H(x) dx
#                   0
#
#       where H(x) is the normalized distribution of image values within mask
    #img = np.array(img)
    #mask = np.array(mask)

    h = img[mask == level]
    print("Extracting SUV% for organ", level)
    
    #histo=[]
    # col=["red", "green", "yellow","pink", "brown", "darkblue", "black", "purple"]
    # i=level
    # plt.hist(h, color= col[i-1],bins=50, range=(0,6))
    # plt.title("Histogram of SUV values for "+ str(level))
    # plt.show()
    suv_x = numpy.percentile(h, p)

    return suv_x

    
def main():  
    
    # Get the current script's directory
    current_dir = os.path.dirname(os.path.abspath(__file__))
    inputCSV = os.path.join(current_dir, "..", 'data', 'cases_retroMManon.csv')                          #input csv file that has the path to each patient image
    
    #Returns concatinated sequnces of evenly spaced numbers over a specified interval - for SUV percentile
    pv=numpy.concatenate((numpy.linspace(10,50,5),
            numpy.linspace(55,80,6),  
            numpy.linspace(82,90,5),
            numpy.linspace(91,100,10)))
    
    pvi=pv.astype(int)      #have to be indigers to be able to append to name
    col_names = ["SUV" + str(x) for x in pvi]
    
    flists = []
                    
    with open(inputCSV, 'r') as inFile:                                         #goes trough the whole csv file and makes an flist (list) 
        cr = csv.DictReader(inFile, lineterminator='\n')
        flists = [row for row in cr]
        
    alldata=[]

    histo=[]
    for idx, entry in enumerate(flists, start=1):                               #
        PET = entry['Image']
        Seg = entry['Mask']
        Visit = entry['Visit']
        ID = entry['ID']
    
        #open PET image and segmentation
        niPET=nibabel.load(PET)
        niSeg=nibabel.load(Seg)
        
        print(" Processing Patient, Visit , (Image: , Mask:)", idx, len(flists), Visit, entry['Image'], entry['Mask'])
        #In nnU-Net (organ segmentation are numbered): 
        #1 liver  #DM
        #2 spleen  #DM
        #3 lungs #DM
        #4 thyroid #DM
        #5 bowel  # DM 16
        #6 pancreas
        #7 bladder
        #8 kidneys
    
        #extract SUV percentiles from PET image with mask - segmentation
        for level in [1,2,3,4,5,6,7,8]:
            v=organ_percentile(niPET.get_fdata(),niSeg.get_fdata(),level,pv)
            mat=[level,ID, Visit]
            mat=numpy.concatenate((mat,v))
            alldata.append(mat)

    # plt.hist(histo, bins=5, range=(0,6))
    # plt.title("Histogram of SUV values for Spleen")
    # plt.show()
            #print("Percentile, value: ", mat)
    #print(alldata)
    # plt.hist(h, bins=50, range=(0,6))
    # plt.show()
    Names=["Organ", "PatientID", "Visit"]
    Names=numpy.concatenate((Names,col_names))
    #make a dataframe to write the data than in excel - that can be read by another script
    df=pd.DataFrame(alldata,columns=Names)
    df.to_excel('..\\data\\percentilesAnon.xlsx', index=None) #path to where you want to save excel
    print("Done.")

if __name__ == '__main__':
  main()