Upload New File

QSM_pipeline/Mix_Mod_CSF.py

+#!/bin/env python
+
+# Description: Script to calculate CSF value for referencing
+#
+# Authors: Chaoyue Wang, Alberto Llera, Benjamin C. Tendler & Karla L. Miller
+#
+# Copyright 2021 University of Oxford
+#
+# Licensed under the Apache License, Version 2.0 (the "License");
+# you may not use this file except in compliance with the License.
+# You may obtain a copy of the License at
+#
+#    http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+#
+
+# This script requires "One_Dim_Mix_Mod" which is available at https://github.com/allera/One_Dim_Mixture_Models
+
+#Add the toolbox to path
+import numpy as np
+import nibabel as nib
+import os
+import sys
+import argparse
+toolbox_path = "./One_Dim_Mix_Mod"
+sys.path.append(os.path.join(os.path.abspath(toolbox_path)))
+
+parser = argparse.ArgumentParser()
+parser.add_argument('-id', dest="ID", type=str, nargs=1)
+argsa = parser.parse_args()
+subID = ' '.join(map(str, argsa.ID))
+
+#load ROI data
+QSM = nib.load('./ROI_QSM_CSF.nii.gz')
+CSF = np.array(QSM.dataobj)
+
+CSF = CSF.ravel()
+data = CSF [CSF != 0]
+
+data_vector=(data - np.mean(data))/np.std(data)
+
+#Define options for the mixture model fit
+Inference ='Method of moments'#'Variational Bayes'#'Method of moments'#'Variational Bayes'  #'Method of moments' OR 'Maximum Likelihood' OR 'Variational Bayes' ML NOT INCLUDED YET
+Number_of_Components=3
+Components_Model=['Gauss','InvGamma','-InvGamma'] #Each component can be Gauss, Gamma, InvGamma, -Gamma, -InvGamma
+init_params=[0,1,5,2,-5,2]
+init_pi=np.ones(3);
+init_pi=np.divide(init_pi,3)
+#init_pi[0]=0.9;init_pi[1]=0.05;init_pi[2]=0.05
+maxits=300
+tol=0.00000001
+opts={'Inference':Inference,'Number_of_Components':Number_of_Components,'Components_Model':Components_Model,
+                                        'init_params':init_params,'maxits':maxits,'tol':tol,'init_pi':init_pi}
+#Define options for the mixture model fit
+
+
+# CALL TO FIT MIXTURE MODEL
+from Mixture_Model_1Dim import Mixture_Model_1Dim     
+Model = Mixture_Model_1Dim(data_vector, opts)
+
+Means=Model['mu1'][0]* np.std(data)+np.mean(data)
+
+f = open("./QSM_REF_CSF_MixMod.txt", "w")
+f.write(str(Means))
+f.close()