VH_Bioinformatics_Shell.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "provenance": []
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "cells": [
    {
      "cell_type": "code",
      "source": [
        "# Necessary pip installs\n",
        "\n",
        "!pip install scanpy\n",
        "!pip install anndata"
      ],
      "metadata": {
        "id": "jMUhvQA9sCaS"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "9Et3ogvVrFyN"
      },
      "outputs": [],
      "source": [
        "# Load all of our necessary imports\n",
        "\n",
        "import scanpy as sc\n",
        "import numpy as np\n",
        "import pandas as pd\n",
        "import matplotlib.pyplot as plt\n",
        "import seaborn as sns\n",
        "import plotly.express as px\n",
        "import plotly.graph_objects as go\n",
        "import anndata as anndata\n",
        "import warnings\n",
        "warnings.filterwarnings('ignore')"
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "# Load the data\n",
        "adata = anndata.read_h5ad('adata_filtered.h5ad')\n",
        "\n",
        "# Source: Aging Fly Cell Atlas (https://hongjielilab.org/afca/)\n",
        "adata"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "FOnyB69ysY1o",
        "outputId": "c5600b2d-b493-46d7-c7af-bceb2332f976"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "AnnData object with n_obs × n_vars = 4619 × 15992\n",
              "    obs: 'tissue', 'sex', 'age', 'sex_age', 'n_genes_by_counts', 'total_counts', 'total_counts_mt', 'pct_counts_mt', 'log1p_n_genes_by_counts', 'log1p_total_counts', 'log1p_total_counts_mt', 'dataset', 'fca_annotation', 'afca_annotation', 'afca_annotation_broad'\n",
              "    uns: 'afca_annotation_colors'\n",
              "    obsm: 'X_pca', 'X_tsne', 'X_umap'"
            ]
          },
          "metadata": {},
          "execution_count": 107
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "# Perform differential gene analysis (FILL OUT)\n",
        "\n",
        "\n",
        "# Extract results into a DataFrame (FILL OUT)\n",
        "result = adata.uns['rank_genes_groups']\n",
        "pvals = None\n",
        "logfc = None\n",
        "\n",
        "# Create DataFrame for plotting (FILL OUT)\n",
        "volcano_df_5to30 = pd.DataFrame({\n",
        "    'Gene': None,\n",
        "    'Log2FoldChange': None,\n",
        "    'P_value': None,\n",
        "    '-log10(P_value)': None\n",
        "})\n",
        "\n",
        "# Sort DataFarme by significance\n",
        "volcano_df_5to30 = volcano_df_5to30.sort_values(by='-log10(P_value)', ascending=False)"
      ],
      "metadata": {
        "id": "lXZ4-OR0tYm4"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "# Volcano plot to see significant differently expressed genes and fold change (FILL OUT)\n",
        "fig = px.scatter(\n",
        "\n",
        ")\n",
        "\n",
        "# Add significance threshold line\n",
        "significance_threshold = -np.log10(0.05)\n",
        "fig.add_trace(\n",
        "    go.Scatter(\n",
        "        x=[min(volcano_df_5to30['Log2FoldChange']), max(volcano_df_5to30['Log2FoldChange'])],\n",
        "        y=[significance_threshold, significance_threshold],\n",
        "        mode=\"lines\",\n",
        "        line=go.scatter.Line(color=\"red\", dash=\"dash\"),\n",
        "        showlegend=False\n",
        "    )\n",
        ")"
      ],
      "metadata": {
        "id": "H_P8TjH5uB30"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "# Heat map to show expression\n",
        "\n",
        "# Format the data correctly: data frame for average gene expression by age (FILL OUT)\n",
        "\n",
        "# Create a seaborn heatmap (FILL OUT)\n",
        "plt.figure(figsize=(6, 7))\n",
        "sns.heatmap()\n",
        "\n",
        "# Formatting\n",
        "plt.xlabel('Age (days)', fontsize=18, fontweight='bold')\n",
        "plt.xticks(fontsize=14)\n",
        "\n",
        "plt.ylabel('Gene', fontsize=18, fontweight='bold')\n",
        "plt.yticks(fontsize=14, rotation = 0)\n",
        "\n",
        "cbar = plt.gca().collections[0].colorbar\n",
        "cbar.ax.tick_params(labelsize=14)\n",
        "cbar.set_label('Mean Expression', size=18, weight='bold', labelpad=20)\n",
        "\n",
        "# Show plot\n",
        "plt.show()"
      ],
      "metadata": {
        "id": "LODIhzEi2GDh"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "# Challenge: Violin plot for Myo81F gene\n",
        "\n",
        "import seaborn as sns\n",
        "import matplotlib.pyplot as plt\n",
        "import numpy as np\n",
        "import pandas as pd\n",
        "\n",
        "# Assuming adata contains your gene expression data\n",
        "gene_name = 'Myo81F'\n",
        "ages = [5, 30, 50, 70]\n",
        "gene_data = {'Age': [], 'Expression': []}\n",
        "\n",
        "# Collect gene expression data for Myo81F across ages in long form\n",
        "for age in ages:\n",
        "    age_filter = adata[(np.array(adata.obs['age'], dtype=int) == age)]\n",
        "    temp = age_filter.X[:, np.where(adata.var_names == gene_name)[0][0]].toarray().flatten()\n",
        "\n",
        "    # Add age and expression to gene_data\n",
        "    gene_data['Age'].extend([age] * len(temp))\n",
        "    gene_data['Expression'].extend(temp)\n",
        "\n",
        "# Convert to DataFrame for plotting\n",
        "df_gene = pd.DataFrame(gene_data)\n",
        "\n",
        "# Create a violin plot\n",
        "\n"
      ],
      "metadata": {
        "id": "ZBH3FHM1LdsQ"
      },
      "execution_count": null,
      "outputs": []
    }
  ]
}