{ "cells": [ { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "# Getting started with pyglotaran\n", "\n", "> Pyglotaran is an open-source modeling framework, written in Python and designed for global and target analysis of time-resolved spectroscopy data. It provides a flexible framework for analyzing complex spectrotemporal datasets, supporting a wide range of kinetic models including sequential, parallel, and target analysis schemes.\n", "\n", "This getting-started notebook exists to get you started with using pyglotaran, hopefully in 20 minutes or less.\n", "\n", "What you are viewing now may be a *static* rendering of an otherwise **interactive** notebook. This guide is the most useful if you either follow along in a *new* notebook, or download the original notebook from the repository.\n", "\n", "
\n", "Click here for more details on notebooks.\n", "\n", "A Python notebook, also known as Jupyter Notebook, is an interactive computational environment, where you can run code, explore data and present your results, all in a single file (with the file extension `.ipynb`).\n", "\n", "There are three main ways to run a Jupyter Notebook:\n", "- `[local]` Using Jupyter Notebook Directly\n", "- `[local|github]` Using VS Code with the Python and Jupyter extensions\n", "- `[cloud]` Using Google Colab \n", "\n", "
\n", "\n", "For the purpose of this guide, it is assumed you already know how to work with notebooks; else there are plenty of tutorials online." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Preface\n", "\n", "If you are going through this guide, you most likely have some dataset burning a hole in your pocket. \n", "\n", "Please rest assured, we'll get to that. But **first**, we'd like to take you through a typically modeling workflow." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Working with data \n", "\n", "Let's start with the premise that we already have some data imported in our notebook.\n", "\n", "To make this guide self-contained, we'll make use of some simulated data from the glotaran.testing package" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# For the purpose of illustrating the workflow, we will just use\n", "# some simulated data from the glotaran.testing package.\n", "from glotaran.testing.simulated_data.sequential_spectral_decay import DATASET as my_dataset\n", "\n", "# ending the cell with the variable name will display the content of the variable\n", "my_dataset" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Like all data in `pyglotaran`, the dataset is an [xarray.Dataset](https://docs.xarray.dev/en/stable/api/dataset.html), for which you can find more information on the [xarray hompage](https://docs.xarray.dev/en/stable).\n", "\n", "From the output cell, we can quickly see that the dataset `my_dataset` has the following properties:\n", "- It has the `Dimensions`: `time` and `spectral` (and they *must* be named like that)\n", "- For *these* data the time coordinate starts at `-1.0` and runs until `19.99`, \n", "- For *these* data the spectral coordinate starts at `600` and runs until `699.4`.\n", "- The dataset (currently) has a single Data variable called `data` with (2100 time x 72 spectral)=151200 datapoints, later more variables may be added to the dataset.\n", "\n", "Towards the end of this notebook you will find out how to read in your data and transform it into an xarray.Dataset, but for now, let's continue." ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "## Plotting raw data\n", "\n", "The pyglotaran framework has built-in functionality to create useful plots. They are part of the `pyglotaran_extras` package, and if you followed the installation instructions, you should already have this installed. \n", "\n", "In there, we have a number of plotting functions we can import and use in our notebook. Let's start with `plot_data_overview`. " ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "from pyglotaran_extras import plot_data_overview\n", "\n", "plot_data_overview(my_dataset);" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Isn't it pretty? We see our 2D intensity map with the time coordinate (labeled Time (ps)) along the x-axis and the spectral coordinate (labeled Wavelength (nm)) along the y-axis.\n", "\n", "Just below that we see our singular value decomposition (SVD) of the data matrix, with the first four (4) singular vectors plotted. \n", "- The left most plot (data. LSV) shows the left singular vectors representing the vector's evolution along the time coordinate.\n", "- The right most plot (data. RSV) shows the right singular vectors reflecting the spectral shapes associated with that evolution. \n", "- The middle plot shows the first 10 singular *values*, on a logarithmic scale. \n", "\n", "From this we can deduce the rank of our data matrix is three (3), which is how many decay components we'll need for our model later.\n", "\n", "
\n", "Q: What? How do you know that? And what's a S.V.D.?\n", "\n", "... A: (click here)\n", "\n", "Very simply put the Singular Value Decomposition (SVD) is a mathematical technique used to decompose a matrix (our data) into a number of left and right singular vectors and their corresponding 'weights' (the singular values). It allows us to quickly visualize the 'rank' of the matrix (which gives us a rough approximation of how many decay components we might need in our model to fit the data). \n", "\n", "More precisely, the SVD decomposed a (data) matrix into three other matrices. It is often used in signal processing and statistics to identify patterns in data. Specifically, SVD decomposes a matrix \\(A\\) into three matrices \\(U\\), \\(Σ\\), and \\(V^T\\) such that \\(A = UΣV^T\\). Here, \\(U\\) and \\(V\\) are orthogonal matrices, and \\(Σ\\) is a diagonal matrix containing the singular values. The left singular vectors (columns of \\(U\\)) represent the time coordinates, while the right singular vectors (columns of \\(V\\)) represent the spectral coordinates. The singular values in \\(Σ\\) provide a measure of the importance of each corresponding singular vector.\n", "\n", "Put simply: if we take the first left singular vector (LSV), multiply it by the first singular value (SV), and then multiply by the first right singular vector (RSV), we obtain the first approximation of the data matrix. Repeating this process for the second and third singular vectors provides an even better approximation, especially since the rank of this particular matrix appears to be 3. By summing the products of all the left and right singular vectors, each weighted by their corresponding singular value, we can reconstruct the original data matrix.\n", "\n", "That's more than you need to know about Singular Value Decompositions at this time. \n", "
" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Starting a project\n", "\n", "Once we have decided that our data is good enough to attempt to model it using `pyglotaran` we can start our adventure.\n", "\n", "To start using `pyglotaran` in your analysis, you only have to import the `Project` class and open a project (folder)." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "from glotaran.project import Project\n", "\n", "my_project = Project.open(\"my_project\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "If the project does not already exist this will create a new project and its folder structure for \n", "you. In our case we had only the `models` + `parameters` folders and the `data` + `results` folder\n", "were created when opening the project." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# This so called shell command allows us to 'list' (ls) the content of the project folder.\n", "%ls my_project" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "### Import the data into your project\n", "\n", "As long as your data can be transformed into a well structured xarray Dataset you can import it directly into your project.\n", "\n", "
\n", "Q: Well structured?\n", "\n", "A: ...\n", "By that we mean that your data is of type `xarray.Dataset` with a `data` Data variable (or `xarray.DataArray`) and has the coordinates `time` and `spectral`.\n", " \n", "This will then save your data as `NetCDF` (`.nc`) file into the data folder inside of your project with the name that you gave it (here `my_project/data/my_data.nc`).\n", "
" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "my_project.import_data(my_dataset, dataset_name=\"my_data\")" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "After importing, our `my_project` is aware of the data that we named `my_data` when importing." ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "### Working with models\n", "\n", "> **⚠️** \n", "> Please note that the *exact* way in which a model is defined may still change slightly in future versions of pyglotaran. But no worries, there will always be a clear procedure to upgrade any existing models you may have created in the meantime.\n", "\n", "After importing our data into the project, to analyse them, we need a `model` (or analysis `scheme`).\n", "\n", "If it does not already exists, create a file called `my_model.yaml` in your projects' `models` folder and fill it with the following content.\n", "\n", "> **📝**\n", "> Don't let this file extension (`.yaml` or `.yml`) scare you, it's just another **plain text file**, which you can open with literally any text editor.\n", "\n", "In our case the file already exists, so we can just *show* you the content, which you can then copy paste." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "my_project.show_model_definition(\"my_model\")\n", "# copy the model definition from the output below" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "⬆️ copy into `models`/`my_model.yaml`\n", "\n", "The above model reads (from bottom to top) as:\n", "- We have a dataset named `my_data`, which we model with a single (kinetic) megacomplex `m1`, with initial_concentration vector `input` and instrument response function (IRF) `irf1`. \n", "- The IRF `irf1` is defined as being of type `gaussian` with its center location at `irf.center` and its width to be `irf.width`\n", "- The megacomplex `m1` is composed of just a single kinetic matrix `k1` of type `decay` (short for exponential decays)\n", "- This k_matrix is composed of just 3 rate constants (`kinetic.1`, `kinetic.2`, `kinetic.3`) describing the sequential kinetic scheme: \"`s1`->`s2`->`s3`->ground\"\n", "- It is sequential because the initial_concentration `input` defines all of the input (1) going to `s1`, and none of it going to `s2` or `s3`.\n", "\n", "You can check your model for problems with the `validate` method." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "my_project.validate(\"my_model\")" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "### Working with parameters\n", "\n", "A model by itself is not sufficient, we also need *starting values* for the parameters we define in the model.\n", "\n", "For this, we use a parameters file. Create a file called `my_parameters.yaml` in your `parameters` folder with the following content.\n", "\n", "Again, in our case the file already exists, so we just *show* the content for you to copy." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "my_project.show_parameters_definition(\"my_parameters\")" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "This reads as: \n", "- There are two `input` parameters `input.1` with (fixed) value 1, and `input.0` with (fixed) value 0.\n", "- There are 3 kinetics rates, with starting values `0.51`, `0.31`, `0.11`.\n", "- There are 2 IRF related parameters `irf.center` with starting value `0.31` and `irf.width` with starting value `0.11`.\n", "\n", "All parameters are implicitly 'free', unless specified with `{ \"vary\": False }` to be fixed." ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "You can use the `validate` method, which can check for missing parameters." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "my_project.validate(\"my_model\", \"my_parameters\")" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "Since not all problems in the model can be detected automatically it is wise to visually inspect the model. \n", "\n", "For this purpose, you can just load the model and inspect its markdown rendered version." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "my_project.load_model(\"my_model\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The same way you should inspect your parameters." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "my_project.load_parameters(\"my_parameters\")" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "## Optimizing - fitting our data\n", "\n", "Now we have all the components in place to start our analysis, optimizing the parameters while minimizing the residual, i.e. fitting our model to our data." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "result = my_project.optimize(\"my_model\", \"my_parameters\")\n", "result" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "What you see from the optimization 'log' is that the optimization took about 5 iterations to converge, going from a cost of 11.2k down to 7.56.\n", "\n", "The Optimization Result table gives us a nice 'statistics' overview and if we click on details we can view the optimized model.\n", "\n", "Since we are satisfied that our fit has converged, we'll proceed to look at the outcome in more detail." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Speaking of outcomes, please note that each time you run an optimization the result will be saved in the project's results folder.\n", "\n", "You may occasionally want to clean these up, especially with larger projects or datasets." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# To view the saved results, look at the content of the project's results folder\n", "%ls \"my_project/results\"" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "One way to look at our results in more detail, is to query the optimized_parameters object in the results.\n", "\n", "This will also print (if it can be computed) the standard error for each *estimated* parameter." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "result.optimized_parameters" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "You can further inspect the `result` by accessing its `data` attribute. In our example it only contains a single `my_data` dataset, but it could contain many datasets in a multi-dataset analysis." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "result.data" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Finally, although knowing your optimization converged, and the optimized parameters seem reasonable is half the battle, the real question is of course ... what does it **LOOK** like. Let's plot!" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Visualize the Result\n", "\n", "The results can be visualized in a similar way as the dataset, using a function from the pyglotaran_extras part of the framework, in this case `plot_overview`. " ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "from pyglotaran_extras import plot_overview\n", "\n", "plot_overview(result);" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In this overview you see the following.\n", "\n", "**first row**:\n", "- `[left]` The concentrations corresponding to the components in the kinetic scheme\n", "- `[middle]` The species associated spectra (SAS), in this case equivalent to the evolution associated spectra (EAS)\n", "- `[right]` The decay associated spectra (DAS)\n", "\n", "**second row**:\n", "- `[left]` The residual matrix, with the IRF (dispersion) curve plotted on top\n", "- `[middle]` The normalized SAS\n", "- `[right]` The normalized DAS\n", "\n", "**third row**: the SVD of the *residual* matrix, with\n", "- `[left]` The first 2 components (black, red) of the left singular vectors (time)\n", "- `[left]` The first 2 components (black, red) of the right singular vectors (spectral)\n", "- `[right]` The first 10 singular values on a logarithmic scale. \n", "\n", "**forth row**: the SVD of the *data* matrix, with\n", "- `[left]` The first 4 components (black, red, blue, green) of the left singular vectors (time)\n", "- `[left]` The first 4 components (black, red, blue, green) of the right singular vectors (spectral)\n", "- `[right]` The first 10 singular values on a logarithmic scale. " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Conclusion\n", "\n", "In this guide we showed how to\n", "- plot your data (once it's in the right format)\n", "- start a project, import the data into it\n", "- work with models (and a model `.yaml` file)\n", "- specify starting values for our parameters (in a parameters `.yaml` file)\n", "- analyze your data\n", "- plot the results\n", "\n", "Welcome to the future of global target analysis, we hope you like it here!" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## And now for something completely different\n", "\n", "Right, right. We were going to talk about _your_ data. That .csv file did end up burning a hole in your pocket, didn't it? 😉\n", "\n", "Well, there is a reason we have standerdized around [xarray](https://docs.xarray.dev/en/stable/)'s [Dataset](https://docs.xarray.dev/en/stable/api/dataset.html)s, we know what to expect.\n", "\n", "We don't know that with `.csv`. 🙈\n", "\n", "- Some people save something as .csv, but it is space or tab delimited. \n", "- Some people put in a header, of a single line, two lines, 4 lines. \n", "- Some people put in a footer, of *many* lines. \n", "- Some people pad their data matrix with a zero, some don't. \n", "- Some people include their spectral and time coordinates, some don't.\n", "- Some do, but call it something completely different. \n", "- Some people use monotonic increasing coordinates (as one should), but some don't 😱.\n", "\n", "In short, you can never be sure with `.csv`." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "All of that being said ... if you had some fairly clean .csv data lying around you wanted to import, you could use the pandas library to read the data and then create an xarray Dataset out of it. \n", "\n", "How, is left as an exercise to the reader, but we left some tips below.\n", "\n", "
\n", "Click to reveal tips!\n", "\n", "Assuming a csv file with the timepoints in the first row, and the spectral coordinate (e.g. wavelengths) in the fist column.\n", "\n", "```csv\n", "-2.0,0.0,2.0,10.0,100.0,1000.0\n", "420,0,10,15,5,1,0\n", "520,1,100,25,7,2,0\n", "620,0,25,15,2,1,0\n", "```\n", "\n", "Then this bit of Python code could read in your data.\n", "\n", "```py\n", "import pandas as pd\n", "import numpy as np\n", "import xarray as xr\n", "from pathlib import Path\n", "filepath = Path(r\"file_that_burned_a_hole_in_your_pocket.csv\")\n", "# Load the CSV file into a pandas DataFrame\n", "df = pd.read_csv(filepath, delimiter=',')\n", "# Convert index and columns to numeric, ignoring any non-numeric values\n", "df.index = pd.to_numeric(df.index, errors='coerce')\n", "df.columns = pd.to_numeric(df.columns, errors='coerce')\n", "# Remove any rows or columns that couldn't be converted to numeric\n", "df = df.loc[df.index.notnull(), df.columns.notnull()]\n", "# Extract the coordinates (assuming they were in the .csv file)\n", "timepoints = np.array(df.columns.values).astype(float)\n", "wavelengths = np.array(df.index.values).astype(float)\n", "dataset = xr.DataArray(\n", " df.values.T,\n", " dims=[\"time\", \"spectral\"],\n", " coords={\"time\": timepoints, \"spectral\": wavelengths},\n", " ).to_dataset(name=\"data\")\n", "dataset\n", "```\n", "\n", "```py\n", "# pro tip, try to plot it to see if it was read in correctly.\n", "plot_data_overview(dataset,linlog=True);\n", "```\n", "\n", "
\n", "\n", "At some point we will add a more robust and battle tested `csv_to_dataset` function to the pyglotaran_extras package, but even then the above continues to hold true." ] } ], "metadata": { "kernelspec": { "display_name": "pyglotaran310", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3" } }, "nbformat": 4, "nbformat_minor": 4 }