Matplotlib के साथ 3D वॉल्यूमेट्रिक डेटा देखना

आपमें से अधिकतर लोग सामान्य कैमरों से ली गई इमेज डेटा (वैज्ञानिक साहित्य में जिन्हें अक्सर “प्राकृतिक छवियाँ” कहा जाता है) से परिचित हैं, साथ ही माइक्रोस्कोप या टेलीस्कोप जैसे विशेष उपकरणों से भी। Python में इमेज के साथ काम करते समय, उन्हें दिखाने का सबसे आम तरीका Matplotlib, जो Python की सबसे लोकप्रिय प्लॉटिंग लाइब्रेरी है, के imshow फ़ंक्शन का उपयोग करना है।

इस ट्यूटोरियल में, हम आपको दिखाएँगे कि इस फ़ंक्शन को 3D वॉल्यूमेट्रिक डेटा दिखाने के लिए कैसे बढ़ाया जाए, जिसे आप इमेजों के एक स्टैक के रूप में सोच सकते हैं। मिलकर, ये एक 3D संरचना का वर्णन करते हैं। उदाहरण के लिए, मैग्नेटिक रेज़ोनेंस इमेजिंग (MRI) और कंप्यूटेड टोमोग्राफी (CT) स्कैन मानव शरीर के भीतर की 3D संरचना को मापते हैं; एक्स-रे माइक्रोटोमोग्राफी काँच या धातु मिश्रधातु जैसे पदार्थों के भीतर की 3D संरचना को मापती है; और लाइट-शीट माइक्रोस्कोप जैविक ऊतकों के भीतर फ्लोरोसेंट कणों को मापते हैं।

हम दिखाएँगे कि MRI डेटासेट कैसे डाउनलोड करें और matplotlib का उपयोग करके उसकी स्लाइसेज़ कैसे प्रदर्शित करें। आप सीखेंगे कि

• nibabel लाइब्रेरी की मदद से अपना डेटा कैसे प्राप्त करें, जो NIfTI फ़ाइल फ़ॉर्मेट के लिए रीडर उपलब्ध कराती है। जो लोग इस चरण को छोड़कर सीधे प्लॉटिंग से शुरू करना चाहते हैं, उनके लिए यह ट्यूटोरियल scikit-image से डेटा प्राप्त करने का तरीका भी दिखाएगा।
• अगला, आप matplotlib के इवेंट हैंडलर API के साथ एक पूर्णत: कार्यात्मक स्लाइस व्यूअर की ओर शुरुआती कदम उठाएँगे, और
• आप अपने नए टूल का उपयोग करने के लिए अपने कुछ कोड को फिर से लिखेंगे।
• अंत में, आप देखेंगे कि काम पूरा होने के बाद आप अपने कार्यस्थल को कैसे साफ कर सकते हैं।

लेकिन पहले कुछ बुनियादी बातें: Matplotlib के imshow से इमेज कैसे प्रदर्शित करें।

यदि आप Jupyter Notebook के साथ काम कर रहे हैं, तो आप इंटरैक्टिव matplotlib मोड सक्षम करके शुरू कर सकते हैं:

%matplotlib notebook

अब आप matplotlib इम्पोर्ट कर सकते हैं और कुछ डेटा प्रदर्शित कर सकते हैं। आप scikit-image लाइब्रेरी के data मॉड्यूल में शामिल कुछ उदाहरण डेटा लोड करेंगे:

import matplotlib.pyplot as plt
from skimage import data

astronaut = data.astronaut()
ihc = data.immunohistochemistry()
hubble = data.hubble_deep_field()

# Initialize the subplot panels side by side
fig, ax = plt.subplots(nrows=1, ncols=3)

# Show an image in each subplot
ax[0].imshow(astronaut)
ax[0].set_title('Natural image')
ax[1].imshow(ihc)
ax[1].set_title('Microscopy image')
ax[2].imshow(hubble)
ax[2].set_title('Telescope image');

नोट: जब आप matplotlib को इंटरैक्टिव नोटबुक मोड में चलाते हैं, तो खुली हुई फ़िगर तब तक सक्रिय रहती है जब तक आप उसे ऊपर-दाएँ कोने में पावर प्रतीक का उपयोग करके निष्क्रिय नहीं कर देते। प्रत्येक प्लॉट से आगे बढ़ने से पहले ऐसा करना सुनिश्चित करें।

इन छवियों को 2-आयामी या 2D छवियाँ कहा जाता है क्योंकि वे 2 आयामों: x और y के साथ व्यवस्थित होती हैं, या NumPy की भाषा में, पंक्तियाँ और स्तंभ या r और c।

कुछ छवियाँ 3D होती हैं, क्योंकि उनमें एक अतिरिक्त गहराई आयाम होता है (z, या प्लेन्स)। इनमें मैग्नेटिक रेज़ोनेंस इमेजिंग (MRI) और सीरियल सेक्शन ट्रांसमिशन इलेक्ट्रॉन माइक्रोस्कोपी (ssTEM) शामिल हैं, जिसमें सैंपल को पतली परतों में काटा जाता है, जैसे सलामी, और प्रत्येक स्लाइस की अलग-अलग इमेजिंग की जाती है।

matplotlib में ऐसी छवियों को देखने के लिए, हमें एक स्लाइस चुननी होती है, और केवल उसी स्लाइस को प्रदर्शित करना होता है। चलिए ऑनलाइन उपलब्ध कुछ मुफ़्त MRI डेटा पर इसे आज़माते हैं।

अंतराल: डेटा प्राप्त करना

हम Buchel और Friston, Cortical Interactions Evaluated with Structural Equation Modelling and fMRI (1997) में वर्णित डेटासेट डाउनलोड करने जा रहे हैं। सबसे पहले, हम डेटा डाउनलोड करने के लिए एक अस्थायी डायरेक्टरी बनाते हैं। हमें अपना विश्लेषण पूरा होने पर इसे हटाना याद रखना चाहिए! यदि आप इस डेटासेट को बाद के लिए रखना चाहते हैं, तो d को अपनी पसंद की किसी स्थायी डायरेक्टरी लोकेशन में बदल दें।

import tempfile

# Create a temporary directory
d = tempfile.mkdtemp()

अब, डेटा डाउनलोड करें:

import os

# Return the tail of the path
os.path.basename('http://google.com/attention.zip')

from urllib.request import urlretrieve

# Define URL
url = 'http://www.fil.ion.ucl.ac.uk/spm/download/data/attention/attention.zip'

# Retrieve the data
fn, info = urlretrieve(url, os.path.join(d, 'attention.zip'))

और इसे zip फ़ाइल से हमारी अस्थायी डायरेक्टरी में निकालें:

import zipfile

# Extract the contents into the temporary directory we created earlier
zipfile.ZipFile(fn).extractall(path=d)

यदि आप फ़ाइल की वास्तविक सामग्री देखें, तो आपको बहुत सारे ‘.hdr’ और ‘.img’ फ़ाइलें मिलेंगी।

# List first 10 files
[f.filename for f in zipfile.ZipFile(fn).filelist[:10]]

ये NIfTI फ़ाइल फ़ॉर्मेट में हैं, और हमें इनके लिए एक रीडर चाहिए। सौभाग्य से, बेहतरीन nibabel लाइब्रेरी ऐसा रीडर उपलब्ध कराती है। सुनिश्चित करें कि आप इसे conda install -c conda-forge nibabel या pip install nibabel से इंस्टॉल करें, और फिर:

import nibabel

अब, अंततः हम अपनी इमेज पढ़ सकते हैं, और देखने के लिए NumPy ऐरे प्राप्त करने के लिए .get_data() मेथड का उपयोग कर सकते हैं:

# Read the image 
struct = nibabel.load(os.path.join(d, 'attention/structural/nsM00587_0002.hdr'))

# Get a plain NumPy array, without all the metadata
struct_arr = struct.get_data()

टिप: यदि आप सीधे MRI डेटा की प्लॉटिंग जारी रखना चाहते हैं, तो निम्नलिखित कोड लाइनों को चलाएँ:

from skimage import io

struct_arr = io.imread("https://s3.amazonaws.com/assets.datacamp.com/blog_assets/attention-mri.tif")

वापस प्लॉटिंग पर

अब उस ऐरे में एक स्लाइस पर नज़र डालते हैं:

plt.imshow(struct_arr[75])

वाह! यह काफी दबा-सा लग रहा है! ऐसा इसलिए है क्योंकि कई MRI में वर्टिकल अक्ष के साथ रिज़ॉल्यूशन हॉरिज़ॉन्टल अक्षों के समान नहीं होता। हम इसे imshow फ़ंक्शन में aspect पैरामीटर पास करके ठीक कर सकते हैं:

plt.imshow(struct_arr[75], aspect=0.5)

लेकिन, चीज़ों को आसान बनाने के लिए, हम डेटा को बस ट्रांसपोज़ कर देंगे और केवल क्षैतिज स्लाइसों को देखेंगे, जिन्हें ऐसे जुगाड़ की ज़रूरत नहीं होती।

struct_arr2 = struct_arr.T
plt.imshow(struct_arr2[34])

सुंदर! बेशक, फिर किसी दूसरी स्लाइस को देखने के लिए, या किसी अलग अक्ष के साथ स्लाइस देखने के लिए, हमें imshow को फिर से कॉल करना होगा:

plt.imshow(struct_arr2[5])

ये सारी कॉलें काफ़ी जल्दी उबाऊ हो जाती हैं। लंबे समय तक, मैं 3D वॉल्यूम्स को Python के बाहर के टूल्स, जैसे ITK-SNAP, से देखता था। लेकिन, जैसा कि पता चलता है, matplotlib व्यूअर में 3D “स्क्रॉलिंग” क्षमताएँ जोड़ना काफ़ी आसान है! इससे हम Python के भीतर ही 3D डेटा का अन्वेषण कर पाते हैं, और डेटा एक्सप्लोरेशन तथा डेटा विश्लेषण के बीच संदर्भ बदलने की ज़रूरत कम हो जाती है।

कुंजी है matplotlib के इवेंट हैंडलर API का उपयोग करना, जो हमें प्लॉट पर किए जाने वाले कार्यों — जिनमें विशेष कुंजी दबाने या माउस बटन क्लिक के जवाब में प्लॉट के डेटा को बदलना भी शामिल है! — को परिभाषित करने देता है।

हमारे मामले में, चलिए की-बोर्ड पर J और K कुंजियों को “पिछली स्लाइस” और “अगली स्लाइस” से बाँधते हैं:

def previous_slice():
    pass

def next_slice():
    pass

def process_key(event):
    if event.key == 'j':
        previous_slice()
    elif event.key == 'k':
        next_slice()

काफ़ी सरल! बेशक, हमें समझना होगा कि इन कार्रवाइयों को वास्तव में कैसे लागू किया जाए और हमें फ़िगर को बताना होगा कि उसे कीबोर्ड प्रेस को प्रोसेस करने के लिए process_key फ़ंक्शन का उपयोग करना चाहिए! दूसरा काम सरल है: हमें फ़िगर कैनवस की mpl_connect मेथड का उपयोग करना है:

fig, ax = plt.subplots()
ax.imshow(struct_arr[..., 43])
fig.canvas.mpl_connect('key_press_event', process_key)

mpl_connect का पूरा प्रलेखन आपको यहाँ मिल सकता है, जिसमें यह भी शामिल है कि आप और किस प्रकार की घटनाओं (जैसे माउस बटन क्लिक) को बाँध सकते हैं।

मुझे ज़्यादा तलाश नहीं करनी पड़ी यह जानने में कि imshow एक AxesImage ऑब्जेक्ट लौटाता है, जो matplotlib के उस Axes ऑब्जेक्ट के “अंदर” रहता है जहाँ सारा ड्रॉइंग होता है, उसकी .images विशेषता में। और यह ऑब्जेक्ट एक सुविधाजनक set_array मेथड देता है जो प्रदर्शित की जा रही इमेज डेटा को बदल देता है! तो, हमें बस ये करना है:

• किसी मनचाहे इंडेक्स को प्लॉट करें, और उस इंडेक्स को सहेजें, शायद रनटाइम के दौरान Axes ऑब्जेक्ट पर एक अतिरिक्त विशेषता के रूप में।
• next_slice और previous_slice फ़ंक्शन उपलब्ध कराएँ जो इंडेक्स बदलें और 3D वॉल्यूम की संबंधित स्लाइस सेट करने के लिए set_array का उपयोग करें।
• नई डेटा के साथ फ़िगर को फिर से ड्रॉ करने के लिए फ़िगर कैनवस की draw मेथड का उपयोग करें।

def multi_slice_viewer(volume):
    fig, ax = plt.subplots()
    ax.volume = volume
    ax.index = volume.shape[0] // 2
    ax.imshow(volume[ax.index])
    fig.canvas.mpl_connect('key_press_event', process_key)

def process_key(event):
    fig = event.canvas.figure
    ax = fig.axes[0]
    if event.key == 'j':
        previous_slice(ax)
    elif event.key == 'k':
        next_slice(ax)
    fig.canvas.draw()

def previous_slice(ax):
    """Go to the previous slice."""
    volume = ax.volume
    ax.index = (ax.index - 1) % volume.shape[0]  # wrap around using %
    ax.images[0].set_array(volume[ax.index])

def next_slice(ax):
    """Go to the next slice."""
    volume = ax.volume
    ax.index = (ax.index + 1) % volume.shape[0]
    ax.images[0].set_array(volume[ax.index])

इसे आज़माएँ!

multi_slice_viewer(struct_arr2)

यह काम कर रहा है! बढ़िया! लेकिन, यदि आप इसे घर पर आज़माएँगे, तो आप नोटिस करेंगे कि K के साथ ऊपर स्क्रॉल करने पर प्लॉट का क्षैतिज स्केल भी दब जाता है। हuh? (यह केवल तब होता है जब आपका माउस इमेज के ऊपर हो।)

हो यह रहा है कि Matplotlib में इवेंट हैंडलर जोड़ने से वे एक-दूसरे के ऊपर चढ़ते जाते हैं। इस मामले में, K x-अक्ष को लॉगरिदमिक स्केल में बदलने का एक बिल्ट-इन कीबोर्ड शॉर्टकट है। यदि हम K का विशेष रूप से उपयोग करना चाहते हैं, तो हमें इसे matplotlib के डिफ़ॉल्ट की मैप्स से हटाना होगा। ये plt.rcParams डिक्शनरी में सूचियों के रूप में रहते हैं, जो matplotlib की सिस्टम-वाइड डिफ़ॉल्ट सेटिंग्स का भंडार है:

plt.rcParams['keymap.<command>'] = ['<key1>', '<key2>']

जहाँ सूची में दी गई किसी भी कुंजी (अर्थात् <key1> या <key2>) को दबाने से <command> निष्पादित हो जाएगा।

इसलिए, हमें एक हेल्पर फ़ंक्शन लिखने की ज़रूरत होगी जो उन कुंजियों को हटा दे जिन्हें हम उपयोग करना चाहते हैं, चाहे वे इस डिक्शनरी में कहीं भी दिखाई दें। (यह फ़ंक्शन अभी matplotlib में मौजूद नहीं है, लेकिन संभवतः यह एक स्वागतयोग्य योगदान होगा!)

def remove_keymap_conflicts(new_keys_set):
    for prop in plt.rcParams:
        if prop.startswith('keymap.'):
            keys = plt.rcParams[prop]
            remove_list = set(keys) & new_keys_set
            for key in remove_list:
                keys.remove(key)

एक पूर्णत: कार्यात्मक स्लाइस व्यूअर

ठीक है, इस नए टूल का उपयोग करने के लिए अपने फ़ंक्शन को फिर से लिखते हैं:

def multi_slice_viewer(volume):
    remove_keymap_conflicts({'j', 'k'})
    fig, ax = plt.subplots()
    ax.volume = volume
    ax.index = volume.shape[0] // 2
    ax.imshow(volume[ax.index])
    fig.canvas.mpl_connect('key_press_event', process_key)

def process_key(event):
    fig = event.canvas.figure
    ax = fig.axes[0]
    if event.key == 'j':
        previous_slice(ax)
    elif event.key == 'k':
        next_slice(ax)
    fig.canvas.draw()

def previous_slice(ax):
    volume = ax.volume
    ax.index = (ax.index - 1) % volume.shape[0]  # wrap around using %
    ax.images[0].set_array(volume[ax.index])

def next_slice(ax):
    volume = ax.volume
    ax.index = (ax.index + 1) % volume.shape[0]
    ax.images[0].set_array(volume[ax.index])

अब, हमें अपने MRI वॉल्यूम में सभी स्लाइसेज़ बिना डिफ़ॉल्ट कीमैप के झंझट भरे हस्तक्षेप के देख पाने चाहिए!

multi_slice_viewer(struct_arr2)

इस विधि की एक अच्छी विशेषता यह है कि यह किसी भी matplotlib बैकएंड पर काम करती है! इसलिए, यदि आप इसे IPython टर्मिनल कंसोल में आज़माते हैं, तो भी आपको ब्राउज़र जैसा ही इंटरैक्शन मिलेगा! और यही बात Qt या Tkinter ऐप पर भी लागू होती है जो matplotlib प्लॉट एम्बेड करता है। यह सरल टूल इसलिए आपको matplotlib की विज़ुअलाइज़ेशन क्षमताओं के इर्द-गिर्द और भी जटिल अनुप्रयोग बनाने देता है।

रुख़्सत होने से पहले

अपने पीछे सफ़ाई करना न भूलें, और अस्थायी डायरेक्टरी (यदि आपने बनाई हो) को हटा दें:

import shutil

# Remove the temporary directory
shutil.rmtree(d)