diff --git a/.idea/.gitignore b/.idea/.gitignore new file mode 100644 index 0000000..13566b8 --- /dev/null +++ b/.idea/.gitignore @@ -0,0 +1,8 @@ +# Default ignored files +/shelf/ +/workspace.xml +# Editor-based HTTP Client requests +/httpRequests/ +# Datasource local storage ignored files +/dataSources/ +/dataSources.local.xml diff --git a/.idea/blog.iml b/.idea/blog.iml new file mode 100644 index 0000000..8b8c395 --- /dev/null +++ b/.idea/blog.iml @@ -0,0 +1,12 @@ + + + + + + + + + + \ No newline at end of file diff --git a/.idea/inspectionProfiles/Project_Default.xml b/.idea/inspectionProfiles/Project_Default.xml new file mode 100644 index 0000000..3dce9c6 --- /dev/null +++ b/.idea/inspectionProfiles/Project_Default.xml @@ -0,0 +1,12 @@ + + + + \ No newline at end of file diff --git a/.idea/inspectionProfiles/profiles_settings.xml b/.idea/inspectionProfiles/profiles_settings.xml new file mode 100644 index 0000000..105ce2d --- /dev/null +++ b/.idea/inspectionProfiles/profiles_settings.xml @@ -0,0 +1,6 @@ + + + + \ No newline at end of file diff --git a/.idea/misc.xml b/.idea/misc.xml new file mode 100644 index 0000000..dc9ea49 --- /dev/null +++ b/.idea/misc.xml @@ -0,0 +1,4 @@ + + + + \ No newline at end of file diff --git a/.idea/modules.xml b/.idea/modules.xml new file mode 100644 index 0000000..e81aaae --- /dev/null +++ b/.idea/modules.xml @@ -0,0 +1,8 @@ + + + + + + + + \ No newline at end of file diff --git a/.idea/vcs.xml b/.idea/vcs.xml new file mode 100644 index 0000000..35eb1dd --- /dev/null +++ b/.idea/vcs.xml @@ -0,0 +1,6 @@ + + + + + + \ No newline at end of file diff --git a/posts/cryo_ET/.ipynb_checkpoints/references-checkpoint.bib b/posts/cryo_ET/.ipynb_checkpoints/references-checkpoint.bib new file mode 100644 index 0000000..212a05e --- /dev/null +++ b/posts/cryo_ET/.ipynb_checkpoints/references-checkpoint.bib @@ -0,0 +1,42 @@ +@article{Shi2023, +author={Shi, Wei +and Cai, Yongfei +and Zhu, Haisun +and Peng, Hanqin +and Voyer, Jewel +and Rits-Volloch, Sophia +and Cao, Hong +and Mayer, Megan L. +and Song, Kangkang +and Xu, Chen +and Lu, Jianming +and Zhang, Jun +and Chen, Bing}, +title={Cryo-EM structure of SARS-CoV-2 postfusion spike in membrane}, +journal={Nature}, +year={2023}, +month={Jul}, +day={01}, +volume={619}, +number={7969}, +pages={403-409}, +abstract={The entry of SARS-CoV-2 into host cells depends on the refolding of the virus-encoded spike protein from a prefusion conformation, which is metastable after cleavage, to a lower-energy stable postfusion conformation1,2. This transition overcomes kinetic barriers for fusion of viral and target cell membranes3,4. Here we report a cryogenic electron microscopy (cryo-EM) structure of the intact postfusion spike in a lipid bilayer that represents the single-membrane product of the fusion reaction. The structure provides structural definition of the functionally critical membrane-interacting segments, including the fusion peptide and transmembrane anchor. The internal fusion peptide forms a hairpin-like wedge that spans almost the entire lipid bilayer and the transmembrane segment wraps around the fusion peptide at the last stage of membrane fusion. These results advance our understanding of the spike protein in a membrane environment and may guide development of intervention strategies.}, +issn={1476-4687}, +doi={10.1038/s41586-023-06273-4}, +url={https://doi.org/10.1038/s41586-023-06273-4} +} + +@article{mepsi2022, +title = {MEPSi: A tool for simulating tomograms of membrane-embedded proteins}, +journal = {Journal of Structural Biology}, +volume = {214}, +number = {4}, +pages = {107921}, +year = {2022}, +issn = {1047-8477}, +doi = {https://doi.org/10.1016/j.jsb.2022.107921}, +url = {https://www.sciencedirect.com/science/article/pii/S1047847722000910}, +author = {Borja {Rodríguez de Francisco} and Armel Bezault and Xiao-Ping Xu and Dorit Hanein and Niels Volkmann}, +keywords = {Simulations, tomographic reconstruction, Image processing, Quality metrics, cryo-EM}, +abstract = {The throughput and fidelity of cryogenic cellular electron tomography (cryo-ET) is constantly increasing through advances in cryogenic electron microscope hardware, direct electron detection devices, and powerful image processing algorithms. However, the need for careful optimization of sample preparations and for access to expensive, high-end equipment, make cryo-ET a costly and time-consuming technique. Generally, only after the last step of the cryo-ET workflow, when reconstructed tomograms are available, it becomes clear whether the chosen imaging parameters were suitable for a specific type of sample in order to answer a specific biological question. Tools for a-priory assessment of the feasibility of samples to answer biological questions and how to optimize imaging parameters to do so would be a major advantage. Here we describe MEPSi (Membrane Embedded Protein Simulator), a simulation tool aimed at rapid and convenient evaluation and optimization of cryo-ET data acquisition parameters for studies of transmembrane proteins in their native environment. We demonstrate the utility of MEPSi by showing how to detangle the influence of different data collection parameters and different orientations in respect to tilt axis and electron beam for two examples: (1) simulated plasma membranes with embedded single-pass transmembrane αIIbβ3 integrin receptors and (2) simulated virus membranes with embedded SARS-CoV-2 spike proteins.} +} \ No newline at end of file diff --git a/posts/vascularNetworks/VascularNetworks.ipynb b/posts/vascularNetworks/VascularNetworks.ipynb index 4fb0f64..c5d54e2 100644 --- a/posts/vascularNetworks/VascularNetworks.ipynb +++ b/posts/vascularNetworks/VascularNetworks.ipynb @@ -9,7 +9,7 @@ "title: \"Vascular Networks\"\n", "author: \"Ali Fele Paranj\"\n", "date: \"November 5 2024\"\n", - "categories: [MATH 612]\n", + "categories: [Graph theory, Vascular Networks]\n", "bibliography: references.bib\n", "---" ] @@ -19,7 +19,7 @@ "id": "637dc720-8a97-46cb-986b-e240a88f0b7d", "metadata": {}, "source": [ - "## Abstract\n", + "# Abstract\n", "I have introduced some basic concepts of micro-circulation and the vascular networks and how they get created (angiogenesis) in health and disease. Then I discuss some angiogenesis models (Anderson-Chaplain as well as BARW) and use the tools of the geomstats to analyze the loopy structure in these networks. I explained the characteristics of the loopy structures in the networks in terms of the parameters of the model. Furthermore, I consider the time evolution of the graphs created by these networks and how the characterization of the loopy structures change through time in these networks." ] }, @@ -28,7 +28,7 @@ "id": "cf9a4a1c-1f0a-45f9-af4f-ae30efcee890", "metadata": {}, "source": [ - "## Introduction" + "# Introduction" ] }, { @@ -146,7 +146,7 @@ "id": "34ab7e20-4aa1-4e5a-8571-535f3e99518d", "metadata": {}, "source": [ - "## Data Pre-Processing\n", + "# Data Pre-Processing\n", "The data generated by the Anderson-Chaplain model or BARW model are typically matrices (i.e. images) (however, there are certain ways to construct a graph representing the final structure during the simulation). To be able to do any form of computation with the simulated data (like finding the loops, etc) it is convenient to turn the images into a graph.\n", "\t\n", "For the outputs of the Anderson-Chaplain model , since the vessels are more crowded at some regions, we first did a dilation followed by a erosion (using OpenCV) to enhance the images for further processing. We then applied a connectivity-8 kernel to the enhanced images to find the nodes, and finally used the processes information to turn the vascular network in the image into a graph in python (using NetworkX library). See figure below for more details. We also applied the same processing pipeline on the data generated by the BARW model which is shown in the following figure" @@ -256,6 +256,359 @@ "\n" ] }, + { + "cell_type": "markdown", + "id": "566595b5-b9b7-4faf-88f9-52871f27c245", + "metadata": {}, + "source": [ + "# Results\n", + "\n", + "We have performed different goemetrical analysis and clustering methods on the loops tructure of the generated graphs by Branching-Annihilating-Random-Walker model. We need to highlight that assigning meaningful biological explanations was not the focus of this study.\n", + "\n", + "In summary, we did the following improvements to the model using the concepts that we learned in this course:\n", + "* Enhanced Loop Detection Algorithms\n", + "* Multiscale Loop Detection WorkFlow\n", + "* Statistical Analysis of Loops\n", + "* Loop Structure Clustering\n", + "* Geometrical Analysis of Loops\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "id": "1801d2eb-9d95-4825-a232-7c70d72068ea", + "metadata": {}, + "source": [ + "## Enhanced Loop Detection Algorithm\n", + "\n", + "Before, we used to generate .png images of the simulation result (see figures above) and then perform image analysis to detect loops. For instance we have convolving the image with 4-connectivity and 8-connectivity matrices to extract the graph structres present in the images. In the new approch, instead, we managed to record the structre of the network in a NetworkX datastructre. This is not easy task to perform without smart useage of the object oriented programming structure for the code. We organized our code into following classes" + ] + }, + { + "cell_type": "markdown", + "id": "b9ed7492-a782-45fe-8d45-16dc96b3e9f5", + "metadata": {}, + "source": [ + "![](images/plantuml.svg)" + ] + }, + { + "cell_type": "markdown", + "id": "b330ccb1-7469-4f12-acba-2601bb21ea65", + "metadata": {}, + "source": [ + "Using this structure, we can record the graph structure of the generated networks as a NetworkX dataframe. Then we can use some of the built-in functions of this library to get the loops (cycles) of the network. However, since the generated networks are large, finding all of the loops (of all scales) is computationally very costly. Instead, we first found a minimal set of cycles in the graph that forms a basis for the cycles space. I.e. we found the loops that can be combined (by symmetric difference) to generate new loops. The following figure shows the basis loops highlighted on the graph." + ] + }, + { + "cell_type": "markdown", + "id": "fe2b6ffb-eefe-4901-ad88-5c3f1f571a11", + "metadata": {}, + "source": [ + "![](images/improvedDetection.png)" + ] + }, + { + "cell_type": "markdown", + "id": "95d57367-bc00-43f5-9f74-3ce17d2a38ce", + "metadata": {}, + "source": [ + "As mentioned above, detected cycles are the basis cycles. The space of all cycles in a graph form a vector space and the basis cycles is a basis for that space. In other words, these cycles are all the cycles necessary to generate all of the cycles in the graph. The addition operation between two cycles is the symmetric difference of their edge set (or XOR of their edges). We can combined the basis cycles to generate higher level (and lower level) structure as shown below." + ] + }, + { + "cell_type": "markdown", + "id": "c214702e-714f-4a42-8c25-0f4ddb3cf4d5", + "metadata": {}, + "source": [ + "![](images/multiScaleCycles.png)" + ] + }, + { + "cell_type": "markdown", + "id": "05d94db0-9d76-4419-81de-7ab8f8c3d5d7", + "metadata": {}, + "source": [ + "We can also extract and scale all of the loops for further analysis. The following figure shows all the loops in the network" + ] + }, + { + "cell_type": "markdown", + "id": "24f85b90-8ce7-41d5-8e9b-a54745b4b59a", + "metadata": {}, + "source": [ + "![](images/grid_basis.png)" + ] + }, + { + "cell_type": "markdown", + "id": "14160fb8-f30a-417f-99fc-eb5a398e769f", + "metadata": {}, + "source": [ + "The following figures shows some of the loop strucgures that we can get by combining the loops above." + ] + }, + { + "cell_type": "markdown", + "id": "92652eef-bc4c-45e6-948f-f4d14b185d68", + "metadata": {}, + "source": [ + "![](images/grid_Generated.png)" + ] + }, + { + "cell_type": "markdown", + "id": "54e8a189-38ce-4e5b-ad3f-43c24bbfdb53", + "metadata": {}, + "source": [ + "## Statistical Analysis of Loops\n", + "\n", + "The mechanism that generated the vascular networks is an stochastic process (Branching process + Simple Radnom Walk process + local interactions (annihilation)). So we need to use statistical notions to make some observations. In the figure below, the histogram of the cycle length is plotted. The interesting observation is the fact that the number of cycles is exponentially distributed (with respect to the Cycle length). The slope of this line (on log-log plot) can reveal some very important facts about the universality class that our model belongs to. Not only this is very interesting and important from theoretical point of view, but also it can have very useful practical applications. For instance, in comparing the simulated network with real vascualr networks, this slope can be one of the components of comparison." + ] + }, + { + "cell_type": "markdown", + "id": "a6eda2cf-aa9c-471d-a309-29508e91b584", + "metadata": {}, + "source": [ + "![](images/CountVsLengthLogLog.png)" + ] + }, + { + "cell_type": "markdown", + "id": "fd0848bb-b725-405d-82bd-2a2a1306a898", + "metadata": {}, + "source": [ + "Furthremore, it is instructive to study the correlation matrix between some of the features of the loop. " + ] + }, + { + "cell_type": "markdown", + "id": "6fb5a5e4-7266-4a6e-9efd-f9e6b64188f4", + "metadata": {}, + "source": [ + "![](images/FeatureCorrelationMatrix.png)" + ] + }, + { + "cell_type": "markdown", + "id": "bacdd9a2-1d12-48e6-ae56-fb4675c74c21", + "metadata": {}, + "source": [ + "## Geometric Shape Analysis: Fréchet and Hausdorff Distances\n", + "\n", + "\n", + "In geometric shape analysis, comparing cycles involves quantifying their similarity based on the spatial arrangement of points in each cycle. Two widely used measures for such comparisons are the **Fréchet Distance** and the **Hausdorff Distance**. These metrics provide different insights into the relationship between cycles, and their results can be visualized as heatmaps of pairwise distances.\n" + ] + }, + { + "cell_type": "markdown", + "id": "3f026db6-2550-4fd3-a19f-e977672b8b81", + "metadata": {}, + "source": [ + "\n", + "### Fréchet Distance\n", + "\n", + "The **Fréchet Distance** between two curves $ A = \\{a(t) \\mid t \\in [0,1]\\} $ and $ B = \\{b(t) \\mid t \\in [0,1]\\} $ is defined as:\n", + "\n", + "$$\n", + "d_F(A, B) = \\inf_{\\alpha, \\beta} \\max_{t \\in [0,1]} \\| a(\\alpha(t)) - b(\\beta(t)) \\|,\n", + "$$\n", + "\n", + "where:\n", + "\n", + "- $ \\alpha(t) $ and $ \\beta(t) $ are continuous, non-decreasing reparameterizations of the curves $ A $ and $ B $.\n", + "- $ \\| \\cdot \\| $ denotes the Euclidean norm.\n", + "- The infimum is taken over all possible parameterizations $ \\alpha $ and $ \\beta $.\n", + "\n", + "\n", + "#### Interpretation of Heatmap\n", + "\n", + "The heatmap for the Fréchet distance shows the pairwise distances between all cycles. Each entry $ (i, j) $ in the heatmap represents $ d_F(C_i, C_j) $, the Fréchet distance between cycle $ C_i $ and cycle $ C_j $. Key insights include:\n", + "\n", + "- **Small Values:** Cycles with low Fréchet distances are geometrically similar in terms of overall shape and trajectory.\n", + "- **Large Values:** Larger distances indicate significant differences in the geometry or shape of the cycles.\n", + "\n", + "The heatmap highlights clusters of similar cycles and outliers with unique geometries." + ] + }, + { + "cell_type": "markdown", + "id": "a2620b65-c498-4a9e-924f-139dd4ad39ac", + "metadata": {}, + "source": [ + "![](images/FrechetDistance.png)" + ] + }, + { + "cell_type": "markdown", + "id": "bf5ff969-3516-453f-89e4-e0de9ac3f711", + "metadata": {}, + "source": [ + "### Hausdorff Distance\n", + "\n", + "The **Hausdorff Distance** between two sets of points $ A $ and $ B $ is defined as:\n", + "\n", + "$$\n", + "d_H(A, B) = \\max \\{ \\sup_{a \\in A} \\inf_{b \\in B} \\| a - b \\|, \\sup_{b \\in B} \\inf_{a \\in A} \\| b - a \\| \\}.\n", + "$$\n", + "\n", + "This can be broken down into:\n", + "\n", + "- $ \\sup_{a \\in A} \\inf_{b \\in B} \\| a - b \\| $: The maximum distance from a point in $ A $ to the closest point in $ B $.\n", + "- $ \\sup_{b \\in B} \\inf_{a \\in A} \\| b - a \\| $: The maximum distance from a point in $ B $ to the closest point in $ A $.\n", + "\n", + "The Hausdorff distance quantifies the greatest deviation between the two sets of points, considering how well one set covers the other.\n" + ] + }, + { + "cell_type": "markdown", + "id": "15a6642d-00df-4135-96c7-030bf2aaa6f2", + "metadata": {}, + "source": [ + "\n", + "### Interpretation of Heatmap\n", + "\n", + "The heatmap for the Hausdorff distance shows pairwise distances between cycles. Each entry $ (i, j) $ represents $ d_H(C_i, C_j) $, the Hausdorff distance between cycle $ C_i $ and cycle $ C_j $. Key insights include:\n", + "\n", + "- **Small Values:** Indicates that the points of one cycle are closely aligned with the points of another cycle.\n", + "- **Large Values:** Reflects that one cycle has points significantly farther away from the other, suggesting geometric dissimilarity.\n", + "\n", + "The heatmap highlights cycles that are well-aligned (small distances) and those that are far apart in terms of shape.\n" + ] + }, + { + "cell_type": "markdown", + "id": "3460a030-f606-4dfe-b013-828020f99d6d", + "metadata": {}, + "source": [ + "![](images/HausdorffDistance.png)" + ] + }, + { + "cell_type": "markdown", + "id": "2aaa0a44-50da-4755-a64a-cffb89f1d5be", + "metadata": {}, + "source": [ + "### Comparison of Metrics\n", + "\n", + "- **Fréchet Distance:** Sensitive to the ordering of points along the curves, making it suitable for comparing trajectories or continuous shapes.\n", + "- **Hausdorff Distance:** Ignores the order of points and focuses on the maximum deviation between sets, making it useful for analyzing shape coverage.\n", + "\n", + "Both metrics complement each other in analyzing the geometric properties of cycles. While the Fréchet distance emphasizes trajectory similarity, the Hausdorff distance focuses on the extent of shape overlap.\n" + ] + }, + { + "cell_type": "markdown", + "id": "61ec85f8-b661-43d3-bd6e-9661b0a75232", + "metadata": {}, + "source": [ + "# Clustering" + ] + }, + { + "cell_type": "markdown", + "id": "1a282cba-6961-4c5b-b18a-85432111ed72", + "metadata": {}, + "source": [ + "#### Features Used for Clustering\n", + "For each cycle in the networks, the following features are extracted:\n", + "\n", + "* Length: Total length of the cycle, calculated as the sum of the lengths of all edges in the cycle.\n", + "* Area: Enclosed area of the cycle, approximated using the Shoelace formula.\n", + "* Compactness: Ratio of the area to the square of the perimeter (area/perimeter$^2$).\n", + "* Centroid: The geometric center of the cycle, defined as the mean of the x- and y-coordinates of all points in the cycle.\n" + ] + }, + { + "cell_type": "markdown", + "id": "fb1e489b-9872-4bee-b98f-892a956fcf77", + "metadata": {}, + "source": [ + "#### Principal Component Analysis (PCA)\n", + " PCA is used for dimensionality reduction by projecting high-dimensional data onto a lower-dimensional subspace while preserving as much variance as possible. The original features (length, area, compactness, centroid) are normalized and reduced to two dimensions using PCA. The first two principal components capture the maximum variance in the data, allowing clustering methods to operate in a simplified feature space.\n", + " \n", + "* Interpretation: The axes in the PCA-reduced space represent linear combinations of the original features. Clusters in this space indicate groups of cycles with similar characteristics.\n" + ] + }, + { + "cell_type": "markdown", + "id": "71f9c08b-c775-43ab-81db-29cb1cf281e5", + "metadata": {}, + "source": [ + "![](images/PCAClustering1.png)" + ] + }, + { + "cell_type": "markdown", + "id": "85b1d077-4a64-42b0-be3e-bea13b4fabda", + "metadata": {}, + "source": [ + "# Feature Transformation and Dimensionality Reduction in Shape Analysis\n", + "\n", + "\n", + "In analyzing geometric properties of cycles in networks, high-dimensional feature spaces often obscure relationships between cycles. This section describes the use of feature engineering and nonlinear dimensionality reduction techniques (**t-SNE** and **UMAP**) to uncover hidden patterns, simplify visualization, and enhance interpretability. The outputs are 2D visualizations that highlight structural and functional groupings within the data. The following composite features were derived:\n", + "\n", + "- **Area-to-Length Ratio:** Computed as $ \\frac{\\text{Area}}{\\text{Length}} $. This feature differentiates between compact cycles (high ratio) and elongated cycles (low ratio).\n", + "- **Normalized Compactness:** Computed as $ \\frac{\\text{Compactness}}{\\text{Length}} $. This feature accounts for the relative efficiency of cycles in enclosing their area.\n" + ] + }, + { + "cell_type": "markdown", + "id": "1014bae2-b485-4e0e-8aed-c234c57696e8", + "metadata": {}, + "source": [ + "\n", + "## Dimensionality Reduction\n", + "\n", + "Nonlinear dimensionality reduction methods project high-dimensional data into a lower-dimensional space while preserving specific structural properties.\n", + "\n", + "### t-SNE (t-Distributed Stochastic Neighbor Embedding)\n", + "\n", + "**t-SNE** minimizes the divergence between probability distributions over pairwise distances in high-dimensional and low-dimensional spaces. It focuses on preserving local structures (relationships between nearby points) and is particularly effective at uncovering clusters. The key parameters are **Perplexity:** Controls the balance between local and global structure (default: 30), and **Output Dimension:** Reduced to 2D for visualization.\n", + "\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "id": "8fc351d3-db9a-4b1e-86f2-0fa02f85254b", + "metadata": {}, + "source": [ + "![Some notes to interpret the plot: that cycles forming tight clusters share strong similarities in features such as length, area, or compactness. Isolated points (outliers) indicate rare or unique geometries. t-SNE emphasizes local structures, making it ideal for detecting smaller, tightly-knit groups.](images/TSNECycles.png)" + ] + }, + { + "cell_type": "markdown", + "id": "3d04f4c8-e32c-4098-a409-d5343b5d95ac", + "metadata": {}, + "source": [ + "\n", + "### UMAP (Uniform Manifold Approximation and Projection)\n", + "\n", + "**UMAP** approximates the high-dimensional data manifold and optimally preserves both local and global structures. It provides more interpretable embeddings with smooth transitions between clusters. The key parameters are **Number of Neighbors:** Defines the size of the local neighborhood considered for embedding (default: 15), and **Output Dimension:** Reduced to 2D for visualization.\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "id": "6078bff6-3ef3-4d07-b86e-9f6499c2c5a6", + "metadata": {}, + "source": [ + "![Some notes to interpret the plot: UMAP preserves both local and global structures, making it suitable for analyzing large-scale patterns. Transitions between clusters indicate gradual changes in feature space, useful for understanding progression or hierarchy in cycle characteristics. Dense clusters suggest strong feature alignment, while sparse areas highlight feature variability.\n", + "](images/UMap.png)" + ] + }, + { + "cell_type": "markdown", + "id": "9c7a178f-ed54-4989-b0b9-4c5e61108426", + "metadata": {}, + "source": [ + "### Conclusion\n", + "\n", + "We used a stochastic process (Branching Annihilating Random Walker) to generate some random networks (that resembes the vascular networks). Then we translated this structure to a networkX data frame for easier processing. We extracted the cycle basis for the cycle space of the graph and using the symmetric difference operation we generated new cycles (of different scales). Then performed different statistical and geometrical analysis on the shape of the loops in the graph. Also we calculated different features for the graph and used dimnsionality reduction methods to see if we can observe any structures (clusters) in low dimension." + ] + }, { "cell_type": "markdown", "id": "99d20ad1-ec5f-452c-806d-850c75bbbf3e", diff --git a/posts/vascularNetworks/VascularNetworks.qmd b/posts/vascularNetworks/VascularNetworks.qmd new file mode 100644 index 0000000..7948129 --- /dev/null +++ b/posts/vascularNetworks/VascularNetworks.qmd @@ -0,0 +1,349 @@ +--- +title: Vascular Networks +author: Ali Fele Paranj +date: November 5 2024 +categories: + - Graph theory + - Vascular Networks +bibliography: references.bib +jupyter: python3 +--- + +# Abstract +I have introduced some basic concepts of micro-circulation and the vascular networks and how they get created (angiogenesis) in health and disease. Then I discuss some angiogenesis models (Anderson-Chaplain as well as BARW) and use the tools of the geomstats to analyze the loopy structure in these networks. I explained the characteristics of the loopy structures in the networks in terms of the parameters of the model. Furthermore, I consider the time evolution of the graphs created by these networks and how the characterization of the loopy structures change through time in these networks. + +# Introduction + +Vascular network in animals is a complex network of vessels to carry the nutrients to and wastes out from the tissue. Main parts of the blood circulation system is determined during the embryo development (like the location and characterization of the main veins and arteries). However, the exact positioning and hierarchy of the capillary networks is not predetermined genetically, but is an adaptive process. A tissue under hypoxia (experiencing low oxygen) will initiate biochemical signals to recruit new vessels. One of the main signaling molecules is the tumor angiogenic factor (TAF) and as the name suggests it is highly expressed by tumor under oxygen stress (hypoxia). The TAF molecules will diffuse in the interstitial space. Upon interaction with pre-existing vessels, they will digest through the basal membrane of the pre-existing vessels and will activate the endothelial cells (that are lining the inner part lumen of all vessels in the body). The endothelial cells activated by TAF will turn into tip cells which will move up in the gradient of TAF, and they will make branches of vessels which will grow as the tip cells are moving. See figure below for details. + + +![Steps of angiogenesis process. The endothelial cells are "activated" by the TAF molecules and turn into tip cells that move towards the source of TAF (up the gradient) and produce new blood vessel segments along their way. This is called tip-induced angiogenesis. This figure is taken from [@Atale2020]](images/Angiogenesis.png){width=60%} + +There are also other factors important in the process of tip cell movement among which is the concentration of fibronectin in the environment. Fibronectin is a protein that acts like a glue and attaches the cells to intracellular matrix, which mediates the cell movement. In the section related to the angiogenesis models we will discuss how this protein is considered in the models. + +## Anderson-Chaplain Model of Angiogenesis + +Anderson-Chaplain model of angiogenesis \cite{Anderson1998} describes the angiogenesis process considering the factors like TAF and fibronectin. This model contains three variables +$\newcommand{\R}{\mathbb{R}}$ +$\newcommand{\abs}[1]{|#1|}$ + +* $n = n(X,t): \Omega \times \R \to \R$: the endothelial-cell density (per unit area). +* $c = c(X,t): \Omega \times \R \to \R$: the tumor angiogenic factor (TAF) concentration (nmol per unit area). +* $f = f(X,t): \Omega \times \R \to \R$: the fibronectin concentration (nmol per unit area). + +and the time evolution is governed by the following system of PDEs + +\begin{align*} + &\frac{\partial n}{\partial t} = D_n\nabla^2 n - \nabla\cdot(\chi n\nabla c) - \nabla\cdot(\rho n \nabla f), \\ + &\frac{\partial c}{\partial t} = -\lambda n c, \\ + &\frac{\partial f}{\partial t} = \omega n - \mu n f, + \end{align*} + +where $D_n$ is a diffusion constant taking the random movement of tip cells into account, $\chi, \rho$ reflects the strength of the chemotaxis of tip cells due to the gradient of TAF, and fibronectin respectively. Furthermore, $\lambda, \mu$ is the rate at which tip cells consume the TAF and fibronectin respectively, and $\omega$ denotes the production of fibronectin by the tip cells. Note that we assume at the start of the angiogenesis process, we have a steady state distribution of fibronectin and TAF and is not diffusing. This assumption is not entirely true and can be enhanced. + +Here in this report, we will be using the discrete and stochastic variation of this model. For more detail see [@Anderson1998]. See figure below for some example outputs of the model. + + +![Some example output of the Anderson-Chaplain model of angiogenesis using the implementation of the model shared by [@Nardini2021]. We have assumed the source of TAF molecules is located at the right edge of the domain, while the pre-existing parent vessels is located at the left edge of the domain. The strength of the chemotaxis and haptotactic (due to fibronectin) signaling is set to be $\chi = 0.4$, and $\rho = 0.4$.](images/Anderson.png){width=60%,fig-align="left",#fig:anderson} + + + +## Branching-Annihilating Random Walker +The Anderson-Chaplain model of angiogenesis is not the only formulation of this phenomena. A popular alternative formulation is using the notion of branching annihilating random walkers for the to explain the branching morphogenesis of vascular networks. A very detailed discussion on this formulation can be found in [@Hannezo2017,@Ucar2021]. This formulation has been also successful to models a vast variety of tip-driven morphogenesis in mammary-glands, prostate, kidney [@Hannezo2017], lymphatic system [@Ucar2023], neural branching [@Ucar2021], and etc. + +The core idea behind this formulation is to assume that the tip cells undergo a branching-annihilating random walk, i.e. they move randomly in the space, turn into pairs randomly (branching), and as they move they produce new cells (stalk) behind their trails, and finally annihilate if they encounter any of the stalk cells. See figure below: + + + +![The network generated by branching-annihilating process, where the tip cells (orange circles) are doing random walk (not necessarily unbiased random walk) and each generate two random walkers at random times (branching). The tip cells make the stalk cells (the blue lines) along their way and the tip cells annihilate when encounter any of the stalk cells.](images/BARW.png){width=50%,fig-align="center",#BARW} + + +# Data Pre-Processing +The data generated by the Anderson-Chaplain model or BARW model are typically matrices (i.e. images) (however, there are certain ways to construct a graph representing the final structure during the simulation). To be able to do any form of computation with the simulated data (like finding the loops, etc) it is convenient to turn the images into a graph. + +For the outputs of the Anderson-Chaplain model , since the vessels are more crowded at some regions, we first did a dilation followed by a erosion (using OpenCV) to enhance the images for further processing. We then applied a connectivity-8 kernel to the enhanced images to find the nodes, and finally used the processes information to turn the vascular network in the image into a graph in python (using NetworkX library). See figure below for more details. We also applied the same processing pipeline on the data generated by the BARW model which is shown in the following figure + + + +![The loop detection processing pipeline performed on the data set generated by the BARW formulation of the branching morphogenesis. Right: the initial raw synthetic data, Left: detected loop through the loop detection pipeline.](images/BARWProc.png){width=50%, fig-align="center"} + + + +![fig-cap: The pre-processing pipeline the synthetic data by Anderson-Chaplain model. Top left: The synthetic data generated by the model. Top right: A series of dilation and erosion applied on the image to enhance the images for network extraction. Bottom left: Applying a connectivity-8 kernel on the image to get the nodes of the graph and converting the processed information to NetworkX graph data structure. Then we apply the cycle finding algorithm in the NetworkX to detect the cycles. Bottom right: The overlay of the detected cycles on the original image.](images/AndersonProc.png){width=50%,fig-align="center"} + +The new data structure that we store the synthetic data (i.e. as graphs in NetworkX) gives us a tremendous flexibility in working and analyzing the shape of each individual loops. See the following figure for further details on this. This data structure will also come in handy when we study the time evolution of the graphs. + + + +![Converting the synthetic data from binary images to graph in NetworkX gives us a great flexibility in working and analyzing (or even evolving) individual loops. For instance, in the figures above, Left and Right sub-figures represent the individual scaled version of the loops of Anderson-Chaplain and BARW models respectively. Note that the shapes that looks more jagged correspond to the larger loops, whereas the loops with dull boundary correspond to the smaller loops.](images/loops.png){width=50%,fig-align="center"} + +**Study Plan**: In this section we will use the notions of Kendall shape space to study the structure of the loops that appear under different parameters of the model. + +## Time Evolution Of Networks + +Vascular networks are not static structure, but rather the evolve in time in response to the changing metabolic demand of the underlying tissue, as well as the metabolic cost of the network itself, and the overall energy required to pump the fluid through the network (See [@Pries2010,Secomb2012,@Pries2014,Ouarne2021] for more discussion). To put this in different words, the role of vascular networks is to deliver nutrients to the tissue and remove the wastes. To do this, it needs to have a space filling configuration with lots of branches. However, due to the Poiseuille law for the flow of fluids in a tube, the power needed to pump the fluid through the tube scales with $r^{-4}$ where $r$ is the radius of the tube. I.e. smaller vessel segments needs a huge power to pump the blood through them. Thus have a massively branched structure is not an optimal solution. On the other hand, the vascular network consists of cells which requires maintenance as well. Thus the optimized vascular network should have a low volume as well. Because of these dynamics in action, in the angiogenesis process first a mesh of new blood vessels form which later evolve to a more ordered and hierarchical structure in a self-organization process. + + + +![Remodeling of vascular network of chick chorioallantoic membrane. Initially (sub-figure 1) a mesh of vascular networks form. Then (sub-figures 2,3,4), through the remodeling dynamics, a more ordered and hierarchical structure emerges. Images are taken from [@Richard2018].](images/remodeling2.png){width=50%,fig-align="center"} + + + +To determine the time evolution of the vascular network we first need to formulate the problem in an appropriate way. First, we represent a given vascular network with a multi-weighted graph $G=(\mathcal{V},\mathcal{E})$ where $V$ is the set of vertices and $E$ is the edge set. We define the pressure $\mathbf{P}$ on the nodes, the flow $ \mathbf{Q} $ on the edges, and let $C_{i,j}, L_{i,j}$ denote the conductivity of an edge, and $L_{i,j}$ denote the length of the same edge. Given the source and sink terms on the nodes $ \mathbf{q} $, the flow in the edges can be determined by + $$\mathcal{L} \mathbf{P} = \mathbf{q},$$ + where $\mathcal{L}$ is the Laplacian matrix of the graph. For more details on this see \autoref{appndix}. Once we know the pressures on the nodes, we can easily calculate the flow through the edges by + $$\bf{Q} = \bf{C} L^{-1} \bf{\Delta} \bf{P}, \tag{2}$$ + where $C$ is a diagonal matrix of the conductance of the edges, $L$ is the diagonal matrix of the length of each edge, $ \Delta $ is the transpose of the incidence matrix, and $ P $ is the pressure on the nodes. $Q$ is the flow of the edges. Once we know the flow in the edges, we can design evolution law to describe the time evolution of the weights of the edges (which by Poiseuille's is a function of the radius of the vessel segment). The evolution law can be derived by defining an energy functional and moving down the gradient of the energy functional to minimize it, or we can take an ad-hoc method and write a mechanistic ODE for time evolution of the conductances. For the energy functional one can write + $$ E(\mathbf{C}) = \frac{1}{2} \sum_{e\in \mathcal{E}}(\frac{Q_e^2}{C_e} + \nu C_e^\gamma), $$ + where $ \mathcal{E} $ is the edge set of the graph, $ Q_e, C_e $ is the flow and conductance of the edge $ e $, and $ \nu,\gamma $ are parameters. The first term in the sum is of the form ``power=current$ \times $potential'' and reflects the power required to pump the flow, and the second term can be shown that reflects the volume of the total network. We can set + $$ \frac{d \mathbf{C}}{dt} = -\nabla E, $$ + which determines the time evolution of the weights in a direction that reduces the total energy. The steady-state solution of this ODE system is precisely the Euler-Lagrange formulation of the least action principle. Alternatively, one can come up with carefully designed ODEs for the time evolution of the conductances that represents certain biological facts. In particular + $$ \frac{d C_e}{dt} = \alpha |Q_e|^{2\sigma} - b C_e + g $$ + proposed by \cite{Ronellenfitsch2016}, and + $$ \frac{d}{dt} \sqrt{C_e} = F(Q_e) - c\sqrt{C_e}, $$ + proposed by \cite{Almeida2022} has been popular choices. See \autoref{fig:GraphEvolution} for more details. It is important to note that in the simulations shown here, the initial network is a toy network. This can be improved by using any of the vascular network generated by any of the angiogenesis models discussed before. + + +![Time evolution of optimal transport network. A triangulation of a 2D domain is considered to be the graph over which we optimize the flow. The sink term is represented by green dot, while the sources are represented by yellow dots. Different sub-figures show the flow network at different time steps towards converging to the optimal configuration.](images/evolution.png){width=50%,fig-align="center"} + + +# Results + +We have performed different goemetrical analysis and clustering methods on the loops tructure of the generated graphs by Branching-Annihilating-Random-Walker model. We need to highlight that assigning meaningful biological explanations was not the focus of this study. + +In summary, we did the following improvements to the model using the concepts that we learned in this course: +* Enhanced Loop Detection Algorithms +* Multiscale Loop Detection WorkFlow +* Statistical Analysis of Loops +* Loop Structure Clustering +* Geometrical Analysis of Loops + + +## Enhanced Loop Detection Algorithm + +Before, we used to generate .png images of the simulation result (see figures above) and then perform image analysis to detect loops. For instance we have convolving the image with 4-connectivity and 8-connectivity matrices to extract the graph structres present in the images. In the new approch, instead, we managed to record the structre of the network in a NetworkX datastructre. This is not easy task to perform without smart useage of the object oriented programming structure for the code. We organized our code into following classes + +![](images/plantuml.svg) + +Using this structure, we can record the graph structure of the generated networks as a NetworkX dataframe. Then we can use some of the built-in functions of this library to get the loops (cycles) of the network. However, since the generated networks are large, finding all of the loops (of all scales) is computationally very costly. Instead, we first found a minimal set of cycles in the graph that forms a basis for the cycles space. I.e. we found the loops that can be combined (by symmetric difference) to generate new loops. The following figure shows the basis loops highlighted on the graph. + +![](images/improvedDetection.png) + +As mentioned above, detected cycles are the basis cycles. The space of all cycles in a graph form a vector space and the basis cycles is a basis for that space. In other words, these cycles are all the cycles necessary to generate all of the cycles in the graph. The addition operation between two cycles is the symmetric difference of their edge set (or XOR of their edges). We can combined the basis cycles to generate higher level (and lower level) structure as shown below. + +![](images/multiScaleCycles.png) + +We can also extract and scale all of the loops for further analysis. The following figure shows all the loops in the network + +![](images/grid_basis.png) + +The following figures shows some of the loop strucgures that we can get by combining the loops above. + +![](images/grid_Generated.png) + +## Statistical Analysis of Loops + +The mechanism that generated the vascular networks is an stochastic process (Branching process + Simple Radnom Walk process + local interactions (annihilation)). So we need to use statistical notions to make some observations. In the figure below, the histogram of the cycle length is plotted. The interesting observation is the fact that the number of cycles is exponentially distributed (with respect to the Cycle length). The slope of this line (on log-log plot) can reveal some very important facts about the universality class that our model belongs to. Not only this is very interesting and important from theoretical point of view, but also it can have very useful practical applications. For instance, in comparing the simulated network with real vascualr networks, this slope can be one of the components of comparison. + +![](images/CountVsLengthLogLog.png) + +Furthremore, it is instructive to study the correlation matrix between some of the features of the loop. + +![](images/FeatureCorrelationMatrix.png) + +## Geometric Shape Analysis: Fréchet and Hausdorff Distances + + +In geometric shape analysis, comparing cycles involves quantifying their similarity based on the spatial arrangement of points in each cycle. Two widely used measures for such comparisons are the **Fréchet Distance** and the **Hausdorff Distance**. These metrics provide different insights into the relationship between cycles, and their results can be visualized as heatmaps of pairwise distances. + + +### Fréchet Distance + +The **Fréchet Distance** between two curves $ A = \{a(t) \mid t \in [0,1]\} $ and $ B = \{b(t) \mid t \in [0,1]\} $ is defined as: + +$$ +d_F(A, B) = \inf_{\alpha, \beta} \max_{t \in [0,1]} \| a(\alpha(t)) - b(\beta(t)) \|, +$$ + +where: + +- $ \alpha(t) $ and $ \beta(t) $ are continuous, non-decreasing reparameterizations of the curves $ A $ and $ B $. +- $ \| \cdot \| $ denotes the Euclidean norm. +- The infimum is taken over all possible parameterizations $ \alpha $ and $ \beta $. + + +#### Interpretation of Heatmap + +The heatmap for the Fréchet distance shows the pairwise distances between all cycles. Each entry $ (i, j) $ in the heatmap represents $ d_F(C_i, C_j) $, the Fréchet distance between cycle $ C_i $ and cycle $ C_j $. Key insights include: + +- **Small Values:** Cycles with low Fréchet distances are geometrically similar in terms of overall shape and trajectory. +- **Large Values:** Larger distances indicate significant differences in the geometry or shape of the cycles. + +The heatmap highlights clusters of similar cycles and outliers with unique geometries. + +![](images/FrechetDistance.png) + +### Hausdorff Distance + +The **Hausdorff Distance** between two sets of points $ A $ and $ B $ is defined as: + +$$ +d_H(A, B) = \max \{ \sup_{a \in A} \inf_{b \in B} \| a - b \|, \sup_{b \in B} \inf_{a \in A} \| b - a \| \}. +$$ + +This can be broken down into: + +- $ \sup_{a \in A} \inf_{b \in B} \| a - b \| $: The maximum distance from a point in $ A $ to the closest point in $ B $. +- $ \sup_{b \in B} \inf_{a \in A} \| b - a \| $: The maximum distance from a point in $ B $ to the closest point in $ A $. + +The Hausdorff distance quantifies the greatest deviation between the two sets of points, considering how well one set covers the other. + + +### Interpretation of Heatmap + +The heatmap for the Hausdorff distance shows pairwise distances between cycles. Each entry $ (i, j) $ represents $ d_H(C_i, C_j) $, the Hausdorff distance between cycle $ C_i $ and cycle $ C_j $. Key insights include: + +- **Small Values:** Indicates that the points of one cycle are closely aligned with the points of another cycle. +- **Large Values:** Reflects that one cycle has points significantly farther away from the other, suggesting geometric dissimilarity. + +The heatmap highlights cycles that are well-aligned (small distances) and those that are far apart in terms of shape. + +![](images/HausdorffDistance.png) + +### Comparison of Metrics + +- **Fréchet Distance:** Sensitive to the ordering of points along the curves, making it suitable for comparing trajectories or continuous shapes. +- **Hausdorff Distance:** Ignores the order of points and focuses on the maximum deviation between sets, making it useful for analyzing shape coverage. + +Both metrics complement each other in analyzing the geometric properties of cycles. While the Fréchet distance emphasizes trajectory similarity, the Hausdorff distance focuses on the extent of shape overlap. + +# Clustering + +#### Features Used for Clustering +For each cycle in the networks, the following features are extracted: + +* Length: Total length of the cycle, calculated as the sum of the lengths of all edges in the cycle. +* Area: Enclosed area of the cycle, approximated using the Shoelace formula. +* Compactness: Ratio of the area to the square of the perimeter (area/perimeter$^2$). +* Centroid: The geometric center of the cycle, defined as the mean of the x- and y-coordinates of all points in the cycle. + +#### Principal Component Analysis (PCA) + PCA is used for dimensionality reduction by projecting high-dimensional data onto a lower-dimensional subspace while preserving as much variance as possible. The original features (length, area, compactness, centroid) are normalized and reduced to two dimensions using PCA. The first two principal components capture the maximum variance in the data, allowing clustering methods to operate in a simplified feature space. + +* Interpretation: The axes in the PCA-reduced space represent linear combinations of the original features. Clusters in this space indicate groups of cycles with similar characteristics. + +![](images/PCAClustering1.png) + +# Feature Transformation and Dimensionality Reduction in Shape Analysis + + +In analyzing geometric properties of cycles in networks, high-dimensional feature spaces often obscure relationships between cycles. This section describes the use of feature engineering and nonlinear dimensionality reduction techniques (**t-SNE** and **UMAP**) to uncover hidden patterns, simplify visualization, and enhance interpretability. The outputs are 2D visualizations that highlight structural and functional groupings within the data. The following composite features were derived: + +- **Area-to-Length Ratio:** Computed as $ \frac{\text{Area}}{\text{Length}} $. This feature differentiates between compact cycles (high ratio) and elongated cycles (low ratio). +- **Normalized Compactness:** Computed as $ \frac{\text{Compactness}}{\text{Length}} $. This feature accounts for the relative efficiency of cycles in enclosing their area. + + +## Dimensionality Reduction + +Nonlinear dimensionality reduction methods project high-dimensional data into a lower-dimensional space while preserving specific structural properties. + +### t-SNE (t-Distributed Stochastic Neighbor Embedding) + +**t-SNE** minimizes the divergence between probability distributions over pairwise distances in high-dimensional and low-dimensional spaces. It focuses on preserving local structures (relationships between nearby points) and is particularly effective at uncovering clusters. The key parameters are **Perplexity:** Controls the balance between local and global structure (default: 30), and **Output Dimension:** Reduced to 2D for visualization. + + + +![Some notes to interpret the plot: that cycles forming tight clusters share strong similarities in features such as length, area, or compactness. Isolated points (outliers) indicate rare or unique geometries. t-SNE emphasizes local structures, making it ideal for detecting smaller, tightly-knit groups.](images/TSNECycles.png) + + +### UMAP (Uniform Manifold Approximation and Projection) + +**UMAP** approximates the high-dimensional data manifold and optimally preserves both local and global structures. It provides more interpretable embeddings with smooth transitions between clusters. The key parameters are **Number of Neighbors:** Defines the size of the local neighborhood considered for embedding (default: 15), and **Output Dimension:** Reduced to 2D for visualization. + + +![Some notes to interpret the plot: UMAP preserves both local and global structures, making it suitable for analyzing large-scale patterns. Transitions between clusters indicate gradual changes in feature space, useful for understanding progression or hierarchy in cycle characteristics. Dense clusters suggest strong feature alignment, while sparse areas highlight feature variability. +](images/UMap.png) + +### Conclusion + +We used a stochastic process (Branching Annihilating Random Walker) to generate some random networks (that resembes the vascular networks). Then we translated this structure to a networkX data frame for easier processing. We extracted the cycle basis for the cycle space of the graph and using the symmetric difference operation we generated new cycles (of different scales). Then performed different statistical and geometrical analysis on the shape of the loops in the graph. Also we calculated different features for the graph and used dimnsionality reduction methods to see if we can observe any structures (clusters) in low dimension. + +## Appendix + +For a graph, the Laplacian matrix contains the information on the in/out flow of stuff into the nodes. +

+ BARW +

+ +Then the Laplacian matrix is given by + $$ D = \begin{pmatrix} + 2 & 0 & 0 & 0 & 0 \\ + 0 & 4 & 0 & 0 & 0 \\ + 0 & 0 & 2 & 0 & 0 \\ + 0 & 0 & 0 & 2 & 0 \\ + 0 & 0 & 0 & 0 & 2 + \end{pmatrix}, $$ + and the adjacency matrix is given by + $$ A = \begin{pmatrix} + 0 & 1 & 1 & 0 & 0 \\ + 1 & 0 & 1 & 1 & 1 \\ + 1 & 1 & 0 & 0 & 0 \\ + 0 & 1 & 0 & 0 & 1 \\ + 0 & 1 & 0 & 1 & 0 + \end{pmatrix}, $$ + and the Laplacian matrix is given by + $$ L = D -A = + \begin{pmatrix} + 2 & -1 & -1 & 0 & 0 \\ + -1 & 4 & -1 & -1 & -1 \\ + -1 & -1 & 2 & 0 & 0 \\ + 0 & -1 & 0 & 2 & -1 \\ + 0 & -1 & 0 & -1 & 2 + \end{pmatrix}. + $$ +It is straight forward to generalize the notion of Laplacian matrix to the weighed graphs, where the degree matrix $ D $, the diagonal entries will be the sum of all weights of the edges connected to that node, and for the adjacency matrix, instead of zeros and ones, we will have the weights of the connections.. + +There is also another way of finding the Laplacian matrix by using the notion of incidence matrix. To do so, we first need to make our graph to be directed. Any combination of the direction on the edges will do the job and will yield in a correct answer. For instance, consider the following directed graph + +For a graph, the Laplacian matrix contains the information on the in/out flow of stuff into the nodes. +

+ BARW +

+ +Its incidence matrix will be + $$ + M = \begin{pmatrix} + -1 & 1 & 0 & 0 & 0 & 0 \\ + 0 & -1 & 1 & -1 & 0 & -1 \\ + 1 & 0 & -1 & 0 & 0 & 0 \\ + 0 & 0 & 0 & 1 & 1 & 0 \\ + 0 & 0 & 0 & 0 & -1 & 1 \\ + \end{pmatrix} + $$ + The Laplacian matrix can be written as + $$ \mathcal{L} = M M^T. $$ + Note that in the case of the weighed graphs, we will have + $$ \mathcal{L} = M W M^T \tag{1}$$ + where $ W $ is a diagonal matrix containing the weights. These computations can be done easily on the NetworkX. + +The incidence matrix is also very useful in calculating the pressure difference between nodes of a particular edge. Let $\Delta = M^T$. Then given the vector $P$ that contains the pressures on the vertices, then the pressure difference on the edges will be given by $\Delta P$, where $\Delta$ is the transpose of the incidence matrix. This comes in handy when we want to calculate the flow of the edges which will be given by +$$ \bf{Q} = \bf{C} L^{-1} \bf{\Delta} \bf{P}, \tag{2} $$ +where $ C $ is a diagonal matrix of the conductance of the edges, $L$ is the diagonal matrix of the ``length'' of each edge, $\Delta$ is the transpose of the incidence matrix, and $P$ is the pressure on the nodes. $Q$ is the flow of the edges. In this particular example we are assuming that the relation between flow and the pressure difference is $Q_e = C_e (p_i - p_j)/L$. But we can have many other choices. + +Knowing the sources and sinks on the nodes, the pressure can be determined by the Kirchhoff law + $$ \mathcal{L} \bf{P} = \bf{q}, $$ +where the vector $ q $ is the sources and the sinks values for each node. This is the same as solving the \textbf{Poisson equation}. This can also be written in terms of the flow, i.e. + $$ \Delta^T \bf{Q} = \bf{q}. $$ + By $ (2) $ we can write + $$ (\bf{\Delta}^T \bf{C}\bf{L}^{-1}\Delta) \bf{P} = \bf{q}. $$ + Since $ \Delta = M^T $, the expression inside the parentheses is clearly Equation (1). + +Similar to the Poisson equation on the graph which is equivalent Kirchhoff's law, we can solve other types of heat and wave equations on the graph as well. The Laplacian matrix play a key role. +$$ \frac{\partial p}{\partial t} = - \mathcal{L} p + q, $$ + for the heat equation, and +$$ \frac{\partial^2 p}{\partial t^2} = -\mathcal{L}p + q, $$ + for the wave equation. + +## References + + diff --git a/posts/vascularNetworks/images/CountVsLengthLogLog.png b/posts/vascularNetworks/images/CountVsLengthLogLog.png new file mode 100644 index 0000000..8353575 Binary files /dev/null and b/posts/vascularNetworks/images/CountVsLengthLogLog.png differ diff --git a/posts/vascularNetworks/images/FeatureCorrelationMatrix.png b/posts/vascularNetworks/images/FeatureCorrelationMatrix.png new file mode 100644 index 0000000..f74c9a4 Binary files /dev/null and b/posts/vascularNetworks/images/FeatureCorrelationMatrix.png differ diff --git a/posts/vascularNetworks/images/FrechetDistance.png b/posts/vascularNetworks/images/FrechetDistance.png new file mode 100644 index 0000000..367dac2 Binary files /dev/null and b/posts/vascularNetworks/images/FrechetDistance.png differ diff --git a/posts/vascularNetworks/images/HausdorffDistance.png b/posts/vascularNetworks/images/HausdorffDistance.png new file mode 100644 index 0000000..9039e63 Binary files /dev/null and b/posts/vascularNetworks/images/HausdorffDistance.png differ diff --git a/posts/vascularNetworks/images/PCAClustering1.png b/posts/vascularNetworks/images/PCAClustering1.png new file mode 100644 index 0000000..8b8f17d Binary files /dev/null and b/posts/vascularNetworks/images/PCAClustering1.png differ diff --git a/posts/vascularNetworks/images/TSNECycles.png b/posts/vascularNetworks/images/TSNECycles.png new file mode 100644 index 0000000..e8d410e Binary files /dev/null and b/posts/vascularNetworks/images/TSNECycles.png differ diff --git a/posts/vascularNetworks/images/UML.png b/posts/vascularNetworks/images/UML.png new file mode 100644 index 0000000..e69de29 diff --git a/posts/vascularNetworks/images/UMap.png b/posts/vascularNetworks/images/UMap.png new file mode 100644 index 0000000..bc12da7 Binary files /dev/null and b/posts/vascularNetworks/images/UMap.png differ diff --git a/posts/vascularNetworks/images/grid_Generated.png b/posts/vascularNetworks/images/grid_Generated.png new file mode 100644 index 0000000..5e70813 Binary files /dev/null and b/posts/vascularNetworks/images/grid_Generated.png differ diff --git a/posts/vascularNetworks/images/grid_basis.png b/posts/vascularNetworks/images/grid_basis.png new file mode 100644 index 0000000..e9e0492 Binary files /dev/null and b/posts/vascularNetworks/images/grid_basis.png differ diff --git a/posts/vascularNetworks/images/improvedDetection.png b/posts/vascularNetworks/images/improvedDetection.png new file mode 100644 index 0000000..daef091 Binary files /dev/null and b/posts/vascularNetworks/images/improvedDetection.png differ diff --git a/posts/vascularNetworks/images/multiScaleCycles.png b/posts/vascularNetworks/images/multiScaleCycles.png new file mode 100644 index 0000000..db8a3c1 Binary files /dev/null and b/posts/vascularNetworks/images/multiScaleCycles.png differ diff --git a/posts/vascularNetworks/images/plantuml.png b/posts/vascularNetworks/images/plantuml.png new file mode 100644 index 0000000..09edccf Binary files /dev/null and b/posts/vascularNetworks/images/plantuml.png differ diff --git a/posts/vascularNetworks/images/plantuml.svg b/posts/vascularNetworks/images/plantuml.svg new file mode 100644 index 0000000..30bb357 --- /dev/null +++ b/posts/vascularNetworks/images/plantuml.svg @@ -0,0 +1,105 @@ +Mainremod: NoneG: NonegammaList: List[float]networkList: List[Network]shapeAnalysis: ShapeAnalysisrun(): voidplot(): voidShapeAnalysisnetworkList: List[Network]cyclesList: Listrun(): voidgetCycles(): voidnewCyclesFromBasisCycles(): voidNetworkductsList: List[Duct]ductsToBeAdded: List[Duct]G: Graphgamma: floatgenerate(): voidupdate(): voidremodel(): voidturnToGraph(): voidDuctnetwork: Networkparent: DucttailCurve: Listtip: Tipupdate(): voidTipposition: Array[float]velocity: Array[float]isAlive: boolcollisionData: tupleupdate(): voidRemodelG: GraphsinkNodes: ListsourceNodes: Listgamma: floatremodel(): voidMovementModelsangleMove(position: Array[float], velocity: Array[float], dt: float): Array[float]cartesianMove(position: Array[float], velocity: Array[float], dt: float): Array[float]DatatumorRad: floattumorPos: Array[float]UtilsrotMat(theta: float): Array[Array[float]]calculate_curve_length(points: List[Array[float]]): floatusesusesuses \ No newline at end of file