conn_plot.html

<!DOCTYPE HTML>
<!--
	Read Only by HTML5 UP
	html5up.net | @ajlkn
	Free for personal and commercial use under the CCA 3.0 license (html5up.net/license)
-->
<html>
	<head>
		<title>CONN toolbox: connectivity matrix</title>
		<meta charset="utf-8" />
		<meta name="description" content="Raphael Vallat, PhD"/>
		<meta name="author" content="Raphael Vallat, PhD" />
		<meta name="keywords" content="Raphael Vallat, raphaelvallat, CONN toolbox, fMRI, functional connectivity, correlation matrix, matlab, python" />
		<meta name="viewport" content="width=device-width, initial-scale=1" />
		<!--[if lte IE 8]><script src="assets/js/ie/html5shiv.js"></script><![endif]-->
		<link rel="stylesheet" href="assets/css/main.css" />
		<link rel="stylesheet" href="assets/css/academicons.css"/>
		<link rel="shortcut icon" type="image/x-icon" href="favicon.ico">
		<link rel="shortcut icon" type="image/x-icon" href="favicon.ico">

		<!-- Add Google Fonts -->
		<link href="https://fonts.googleapis.com/css?family=Noto+Sans|Noto+Serif|Roboto+Mono" rel="stylesheet">
		<!-- Add Highlight JS -->
		<script src="assets/js/highlight.pack.js"></script>
		<script>hljs.initHighlightingOnLoad();</script>
		<link rel="stylesheet" href="assets/css/highlightjs/atom-one-light.css">
		<!--[if lte IE 8]><link rel="stylesheet" href="assets/css/ie8.css" /><![endif]-->
		<!-- Global site tag (gtag.js) - Google Analytics -->
		<script async src="https://www.googletagmanager.com/gtag/js?id=UA-118198513-1"></script>
		<script>
			window.dataLayer = window.dataLayer || [];
			function gtag() {
				dataLayer.push(arguments);
			}
			gtag('js', new Date());
			gtag('config', 'UA-118198513-1');
		</script>
	</head>
	<body>

		<!-- Header -->
			<section id="header">
				<header>
					<!-- <span class="image avatar"><img src="images/avatar.jpg" alt="" /></span> -->
					<h1 id="logo"><a href="#">Raphael Vallat, PhD</a></h1>
					<p>Senior ML data scientist<br>Oura Ring</p>
				</header>
				<nav id="nav">
					<ul>
						<li><a href="index.html#one">About</a></li>
		        <li><a href="index.html#two">Academics</a></li>
		        <li><a href="index.html#three">Publications</a></li>
		        <li><a href="index.html#four">Media</a></li>
		        <li><a href="index.html#five">Softwares</a></li>
		        <li><a href="index.html#six" class="active">Tutorials</a></li>
		      </ul>
				</nav>
				<footer>
					<ul class="icons">
						<li><a href="https://github.com/raphaelvallat" class="icon fa-github" target="_blank"><span class="label">Github</span></a></li>
						<li><a href="mailto:raphaelvallat9@gmail.com" class="icon fa-envelope"><span class="label">Email</span></a></li>
						<li><a href="https://scholar.google.fr/citations?user=XH8IG2UAAAAJ" style="border:none" class="ai ai-google-scholar ai-1x" target="_blank"></a></li>
						<li><a href="https://twitter.com/RaphaelVallat" class="icon fa-twitter"><span class="label">Twitter</span></a></li>
						<li><a href="https://www.linkedin.com/in/raphaelvallat9/" class="icon fa-linkedin"><span class="label">Twitter</span></a></li>
					  </ul>
				</footer>
			</section>

		<!-- Wrapper -->
			<div id="wrapper">

				<!-- Main -->
					<div id="main">

					<div class="container">


						<header>
							<h4>CONN Toolbox</h4>
							<h3>Extract and plot ROI-to-ROI connectivity matrix</h3>
							<p>February 2017</p>
						</header>


						<p align="justify">The <a href="https://www.nitrc.org/projects/conn" target="_blank"> CONN toolbox </a> is a very efficient MATLAB-based tool to perform functional connectivity analysis and I've been an enthusiastic user ever since I discovered it. It provides a robust analysis pipeline from raw files to second-level adavanced statistical results as well as a very nice and intuitive GUI - making it one of the best fcMRI toolbox for both beginner and expert. I worked a great deal with it for ROI-to-ROI analysis of within and between network connectivity. One feature that is not (yet) implemented in CONN is the possibility to plot a graphical representation of network connectivity matrix. The purpose of this tutorial is therefore to show how to extract (and export to CSV) connectivity matrix from CONN results folder and then plot them, using either Matlab or Python.

						<ul class="feature-icons"><li class="fa-code">A full version of the MATLAB / Python scripts used can be found on my <a href="https://github.com/raphaelvallat/conn_tools" target="_blank">GitHub repository</a></li></ul>


						<span class="image fit"><img src="images/tutorials/conn_correl_mat/correl_mat.png" alt="" />
						<figcaption>Salience network connectivity matrices (beta, T and p-values)</figcaption>
						</span>

						<hr />

						<h4>ROI.mat structure</h4>

						<p align="justify">CONN ROI.mat file is created when clicking "Results Explorer" in the bottom left of the CONN 2nd-level interface. This MAT structure can be found in the folder <i>/results/secondlevel/FIRSTLEVEL_NAME/GROUP/SESSION/</i> and contains the following fields <a href="https://www.nitrc.org/forum/message.php?msg_id=15500" target="_blank">(see this post)</a></p>

						<li><b>names:</b> names of the source ROIs</li>
						<li><b>h:</b> beta values displayed in CONN 2nd-level interface, which corresponds to the  average Fischer transformed pairwise correlations of specified contrast</li>
						<li><b>F:</b> statistical values, depending on the statistical test (see statsname field)</li>
						<li><b>p:</b> one-tailed p-values</li>

						<hr />

						<h3>1) Using MATLAB </h3>
						<h4>Extract connectivity matrix</h4>
						<p>The following code extract all the pairwise correlations in the ROI.mat file </p>

<pre><code class="matlab">% Define analysis path
wdir 		= 'C:/Users/Raphael/Desktop/These/conn_example/results/secondlevel';
corr_net 	= 'Salience';
corr_group 	= 'AllSubjects';
corr_run 	= 'rest';
corr_folder = [ wdir '/' corr_net '/' corr_group '/' corr_run '/' ];

load([corr_folder 'ROI.mat']);

numROI  	= size(ROI, 2);
corr_name 	= ROI(1).names(1:numROI);

corr_h   	= [];   % Beta value
corr_F  	= [];   % Statistic value
corr_p   	= [];   % One-tailed p value

for i = 1:numROI
	corr_h = [ corr_h ; ROI(i).h(1:numROI) ];
	corr_F = [ corr_F ; ROI(i).F(1:numROI) ];
	corr_p = [ corr_p ; ROI(i).p(1:numROI) ];
end
</code></pre>

						<hr />

						<h4>Export to CSV</h4>
						<p>Now we can export it to CSV file (by default in the same folder as ROI.mat file</p>
<pre><code class="matlab">% array2table does not work with '-' in var names
corr_name2 	= strrep(corr_name, '-', '_');
T_h        	= array2table(corr_h, 'RowNames', corr_name2, 'VariableNames', corr_name2);
T_F        	= array2table(corr_F, 'RowNames', corr_name2, 'VariableNames', corr_name2);
T_p        	= array2table(corr_p, 'RowNames', corr_name2, 'VariableNames', corr_name2);

writetable( T_h, [corr_folder 'beta_' corr_net '_' corr_group '_' corr_run '.csv'], 'WriteVariableNames', true, 'WriteRowNames', true, 'delimiter', 'semi' );
writetable( T_F, [corr_folder 'F_' corr_net '_' corr_group '_' corr_run '.csv'], 'WriteVariableNames', true, 'WriteRowNames', true, 'delimiter', 'semi');
writetable( T_p, [corr_folder 'p_' corr_net '_' corr_group '_' corr_run '.csv'], 'WriteVariableNames', true, 'WriteRowNames', true, 'delimiter', 'semi');
</code></pre>

						<hr />

						<h4>Statistics</h4>
						<p>The following lines compute uncorrected, bonferroni and FDR-corrected p-values</p>
<pre><code class="matlab">anti_corr_net = False;

% Case two or several anti-correlated networks (ex: Default and Dorsal Attention)
if anti_corr_net
    tail      = 'two-sided';
    corr_p    = 2*min(corr_p, 1-corr_p);

% Case single network
else
    tail      = 'one-sided';

end

% Compute Bonferroni and FDR correction
% conn_fdr function is in conn main folder
alpha_bonf            		= 0.05 / ((numROI)*(numROI-1)/2);
vector_fdr                      = nonzeros(triu(corr_p)');
vector_fdr(isnan(vector_fdr))   = [];
corr_p_fdr            		= conn_fdr(vector_fdr);

% WRITE OUTPUT
fprintf('\nANALYSIS INFO');
fprintf('\n--------------------------------------');
fprintf(['\nNetwork:\t ' corr_net]);
fprintf(['\nGroup:\t\t ' corr_group]);
fprintf(['\nRun:\t\t ' corr_run]);
fprintf('\nSTATISTICS');
fprintf('\n--------------------------------------');
fprintf([ '\n' num2str(numROI) ' x ' num2str(numROI-1) ' ROIs matrix ; ' tail]);
fprintf([ '\np-uncorrected:\t\t\t\t\t ' num2str(numel(corr_p(corr_p <= 0.05))/2)  ]);
fprintf([ '\np-bonferroni (alpha = ' num2str(round(alpha_bonf, 5)) '):\t ' num2str(numel(corr_p(corr_p <= alpha_bonf))/2) ]);
fprintf([ '\np-FDR corrected:\t\t\t\t ' num2str(numel(corr_p_fdr(corr_p_fdr <= 0.05))) ]);
fprintf('\n--------------------------------------\n');
</code></pre>

<h5>MATLAB Output:</h5>
<pre>
ANALYSIS INFO
--------------------------------------
Network:	 DMN
Group:		 Controls(1).Patients(-1)
Run:		 rest
STATISTICS
--------------------------------------
5 x 4 ROIs matrix ; one-sided
p-uncorrected:                  4
p-bonferroni (alpha = 0.005):   1
p-FDR corrected:                1
--------------------------------------
</pre>
						<hr />
						<h4>Plotting brain connectivity matrix</h4>
						<p align="justify">Once we extracted and exported the connectivity matrix in text format, we can plot it using Matlab imagesc command. I will detail here for beta values (corr_h matrix) but if you want to plot F or p values you can check the full script on my <a href="https://github.com/raphaelvallat/conn_tools" target="_blank">GitHub repository</a></p>

<pre><code class="matlab">% Remove upper triangle + diagonal
corr = tril(corr, -1)

% Start plotting
fig = figure;
set(gcf,'Units','inches', 'Position',[0 0 6 4])

im = imagesc(corr, clim );

clim = [ 0 1 ];
colormap(flipud(hot(10)));

h = colorbar('eastoutside');
xlabel(h, 'h', 'FontSize', 14);


% Title, axis
title('Salience Network', 'FontSize', 14);
set(gca, 'XTick', (1:numROI));
set(gca, 'YTick', (1:numROI));
set(gca, 'Ticklength', [0 0])
grid off
box off

% Labels
set(gca, 'XTickLabel', corr_name, 'XTickLabelRotation', 0);
set(gca, 'YTickLabel', corr_name);
</code></pre>

						<span class="image fit"><img src="images/tutorials/conn_correl_mat/Salience.png" alt="" /></span>

						<hr />

						<p>With some adaptation, the previous code also works for between-network connectivity matrix  (blue denotes negative correlations)</p>

						<span class="image fit"><img src="images/tutorials/conn_correl_mat/ACN.png" alt="" /></span>

						<hr />

						<h3>2) Using Python </h3>
						<p>The previous connectivity matrices can also be extracted and plot using Python (requires <a href="http://www.numpy.org/" target="_blank">numpy</a>, <a href="http://scipy.org/" target="_blank">scipy</a>, <a href="http://pandas.pydata.org/" target="_blank">pandas</a> and <a href="http://seaborn.pydata.org/index.html" target="_blank">seaborn</a> packages). Apart form being free, Python has several advantages compared to Matlab, with notably the possibility to emphasize networks separation (see example below), and to annotate values in each cell of the heatmap plot. The only trick is in the importation and conversion of the original MAT file to a Pandas dataframe. </p>

						<span class="image fit"><img src="images/tutorials/conn_correl_mat/python_corr_mat.png" alt="" /></span>
<h5>Code:</h5>
<pre><code class="python">import numpy as np
import pandas as pd
import seaborn as sns
import matplotlib.pyplot as plt
import os
import scipy.io

sns.set(style="white", context='paper', font_scale=1, font='monospace')

def extract_conn_correl_mat():

    # Analysis Information
    wdir             = 'C:/Users/Raphael/Desktop/These/CONN_Club_Neuro/Conn_ClubNeuro_Example/results/secondlevel'
    correl_net       = 'Salience'
    correl_group     = 'AllSubjects'
    correl_run       = 'rest'
    correl_type      = 'p'
    correl_folder    = os.path.join(wdir, correl_net, correl_group, correl_run )

    # Load .mat file
    mdata           = scipy.io.loadmat(os.path.join(correl_folder, 'ROI.mat'))
    mdata           = mdata['ROI']
    mdtype          = mdata.dtype
    ndata           = {n: mdata[n] for n in mdtype.names}
    mcol            = ['names', 'h', 'F', 'p']

    mcorr           = pd.DataFrame(np.concatenate([ndata[c] for c in mcol], axis=0)).transpose()
    mcorr.columns   = mcol
    mcorr.names     = np.concatenate(mcorr.names[0][0])
    mcorr.names     = mcorr.names.str.split('.').str.get(0)

    # Extract selected values
    if correl_type == 'beta':
        corr = pd.DataFrame(np.concatenate([mcorr.h[c] for c in mcorr.h.keys()], axis=0 ), index=mcorr.names)
    elif correl_type == 'F':
        corr = pd.DataFrame(np.concatenate([mcorr.F[c] for c in mcorr.F.keys()], axis=0 ), index=mcorr.names)
    elif correl_type == 'p':
        corr = pd.DataFrame(np.concatenate([mcorr.p[c] for c in mcorr.p.keys()], axis=0 ), index=mcorr.names)

    corr            = corr.iloc[:, 0:corr.shape[0]]
    corr.columns    = mcorr.names

    # Export to csv
    corr.to_csv(os.path.join(correl_folder, correl_type + '_' + correl_net + '_' + correl_group + '_' + correl_run + '_corr_mat.csv'), sep=';', decimal='.')

    # Plot
    plot_correl_matrix(corr, correl_type, correl_net, correl_group, correl_run, correl_folder )

def plot_correl_matrix(corr, correl_type, correl_net, correl_group, correl_run, correl_folder ):

    # Mask diagonal
    mask = np.zeros_like(corr, dtype=np.bool)
    mask[np.triu_indices_from(mask, k=1)] = True

    # Define plot properties
    if correl_type == 'F':
        vmin    = 0
        vmax    = 10
        annot   = False
        cmap    = "YlOrRd"

    elif correl_type == 'beta' :
        vmin    = 0
        vmax    = 1
        annot   = True
        cmap    = "YlOrRd"

    elif correl_type == 'p' :
        vmin    = 0
        vmax    = 0.1
        annot   = False
        cmap    = "YlOrRd_r"

    f, ax = plt.subplots()
    sns.heatmap(corr, mask=mask, vmin=vmin,  vmax=vmax, square=True, cmap=cmap, annot=annot, cbar=True, xticklabels=False, yticklabels=True, linewidths=.0 )

    plt.xticks(rotation=0)
    plt.ylabel('')
    plt.xlabel('')
    plt.title(correl_type + '_' + correl_net + '_' + correl_group + '_' + correl_run)

	# Emphasize networks separation
	#    networks = np.array(pd.read_csv(os.path.join(correl_folder, 'networks_level_values.csv')))
	#    for i, network in enumerate(networks):
	#        if i and network != networks[i - 1 ]:
	#            ax.axhline(len(networks) - i, c="w")
	#            ax.axvline(i, c="w")
	#            f.tight_layout()</font>

    plt.savefig(os.path.join(correl_folder,  correl_type + '_' + correl_net + '_' + correl_group + '_' + correl_run + '.png'), dpi=300)

if __name__ =='__main__':

    extract_conn_correl_mat()
</code></pre>

				</div>

				<!-- Footer -->
					<section id="footer">
						<div class="container">
							<ul class="copyright">
								<li>&copy; Raphael Vallat </li><li>Design: <a href="http://html5up.net">HTML5 UP</a></li><li>Banner image: Nicholas Roerich</li>
							</ul>
						</div>
					</section>

			</div>

		<!-- Scripts -->
			<script src="assets/js/jquery.min.js"></script>
			<script src="assets/js/jquery.scrollzer.min.js"></script>
			<script src="assets/js/jquery.scrolly.min.js"></script>
			<script src="assets/js/skel.min.js"></script>
			<script src="assets/js/util.js"></script>
			<!--[if lte IE 8]><script src="assets/js/ie/respond.min.js"></script><![endif]-->
			<script src="assets/js/main.js"></script>

	</body>
</html>