{
"cells": [
{
"cell_type": "markdown",
"id": "078900c6-d71d-44f0-9e76-cc5ce66f0ca9",
"metadata": {},
"source": [
"# Logistic regression: fit\n",
"\n",
"Exploring different fitting techniques and using MLE asymptotic approx of parameter distribution to estimate model CI intervals."
]
},
{
"cell_type": "markdown",
"id": "e18cec5e",
"metadata": {},
"source": [
"## Common"
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "3b351af7",
"metadata": {},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
"import numpy as np\n",
"import scipy\n",
"\n",
"import pandas as pd\n",
"import seaborn as sns\n",
"import os\n",
"\n",
"from sklearn import linear_model\n",
"\n",
"# our libs\n",
"import data_utils\n",
"import logit_utils"
]
},
{
"cell_type": "markdown",
"id": "e0de7f02-942b-42ff-874c-3c0e2ac1e0dc",
"metadata": {},
"source": [
"## Data"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "3373f87f",
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"results_path = \"./results\"\n",
"if not os.path.exists(results_path): os.makedirs(results_path)\n",
"\n",
"data_path = \"../../data/\"\n",
"suv_filename = os.path.join(data_path, \"suv_percentilesSLOthenUWM.mat\")\n",
"flags_filename = os.path.join(data_path,\"flags_combined.mat\")\n",
"normal_range_filename = os.path.join(data_path, \"normal_range.mat\")"
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "811d990d",
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"\"\"\"\n",
"Loading all data\n",
"\"\"\"\n",
"suv_dict = scipy.io.loadmat(suv_filename)\n",
"flags_dict = scipy.io.loadmat(flags_filename)\n",
"normal_range_dict = scipy.io.loadmat(normal_range_filename)"
]
},
{
"cell_type": "code",
"execution_count": 4,
"id": "4a292c31",
"metadata": {},
"outputs": [],
"source": [
"# get data\n",
"perc = 95\n",
"organ = \"lung\"\n",
"\n",
"x0, y = data_utils.get_data(organ, perc, suv_dict, flags_dict)\n",
"\n",
"labs = [\"plain\", \"log\"]\n",
"xs = [x0, np.log(x0)]"
]
},
{
"cell_type": "markdown",
"id": "35e9bb7b-5a50-4cae-9af4-261c11453dd4",
"metadata": {},
"source": [
"## Fit\n",
"\n",
"Ref: logistic regression documentation\n",
"* https://scikit-learn.org/stable/modules/linear_model.html#logistic-regression\n",
"* https://stats.stackexchange.com/questions/186830/what-is-scikit-learns-logisticregression-minimizing"
]
},
{
"cell_type": "code",
"execution_count": 5,
"id": "239ae789",
"metadata": {
"tags": []
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"parss = array([[-11.68319035, 5.68708926],\n",
" [ -7.52531005, 10.73177719]])\n"
]
}
],
"source": [
"# fits\n",
"logr = linear_model.LogisticRegression(penalty = None)\n",
"parss = np.array([logit_utils.logit_poly_fit(logr, x, y, degree = 1) for x in xs])\n",
"print(f\"{parss = }\")"
]
},
{
"cell_type": "code",
"execution_count": 6,
"id": "9cf6fec1",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"cov_parss = array([[[12.31854544, -6.71118734],\n",
" [-6.71118734, 3.81820749]],\n",
"\n",
" [[ 5.04803783, -7.92922372],\n",
" [-7.92922372, 13.99363052]]])\n"
]
}
],
"source": [
"# asymptotic covariance matrix of parameters\n",
"cov_parss = np.array([logit_utils.logit_poly_cov(x, pars) for x, pars in zip(xs, parss)])\n",
"print(f\"{cov_parss = }\")"
]
},
{
"cell_type": "code",
"execution_count": 7,
"id": "2359d3c1",
"metadata": {},
"outputs": [],
"source": [
"# defining two-sided confidence intervals\n",
"alpha = 0.05 # significance level\n",
"probs = [alpha/2, 1 - alpha/2]"
]
},
{
"cell_type": "code",
"execution_count": 8,
"id": "4ba4d5d4",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"pars_CIs = array([[[-18.5622299 , 1.85727375],\n",
" [ -4.80415081, 9.51690477]],\n",
"\n",
" [[-11.92892555, 3.3999319 ],\n",
" [ -3.12169455, 18.06362248]]])\n"
]
}
],
"source": [
"# CI of params (assuming asymptotic distr of parameters)\n",
"pars_CIs = np.array([logit_utils.logit_poly_pars_quantiles_normal(probs, pars, cov_pars) for pars, cov_pars in zip(parss, cov_parss)])\n",
"print(f\"{pars_CIs = }\")"
]
},
{
"cell_type": "code",
"execution_count": 9,
"id": "96acc3ab",
"metadata": {},
"outputs": [
{
"data": {
"application/vnd.microsoft.datawrangler.viewer.v0+json": {
"columns": [
{
"name": "('scale', 'coef')",
"rawType": "object",
"type": "unknown"
},
{
"name": "value",
"rawType": "float64",
"type": "float"
},
{
"name": "LCL",
"rawType": "float64",
"type": "float"
},
{
"name": "UCL",
"rawType": "float64",
"type": "float"
}
],
"ref": "f3aae3c6-f11a-4234-ad63-287af6c89f82",
"rows": [
[
"('plain', 'b0')",
"-11.683190353932442",
"-18.562229897321515",
"-4.804150810543373"
],
[
"('plain', 'b1')",
"5.687089258005262",
"1.8572737464039784",
"9.516904769606544"
],
[
"('log', 'b0')",
"-7.5253100481979",
"-11.928925546652636",
"-3.1216945497431663"
],
[
"('log', 'b1')",
"10.731777187587987",
"3.399931895837665",
"18.063622479338306"
]
],
"shape": {
"columns": 3,
"rows": 4
}
},
"text/html": [
"\n",
"\n",
"
\n",
" \n",
" \n",
" | \n",
" | \n",
" value | \n",
" LCL | \n",
" UCL | \n",
"
\n",
" \n",
" | scale | \n",
" coef | \n",
" | \n",
" | \n",
" | \n",
"
\n",
" \n",
" \n",
" \n",
" | plain | \n",
" b0 | \n",
" -11.683190 | \n",
" -18.562230 | \n",
" -4.804151 | \n",
"
\n",
" \n",
" | b1 | \n",
" 5.687089 | \n",
" 1.857274 | \n",
" 9.516905 | \n",
"
\n",
" \n",
" | log | \n",
" b0 | \n",
" -7.525310 | \n",
" -11.928926 | \n",
" -3.121695 | \n",
"
\n",
" \n",
" | b1 | \n",
" 10.731777 | \n",
" 3.399932 | \n",
" 18.063622 | \n",
"
\n",
" \n",
"
\n",
"
\n",
"\n",
"
\n",
" \n",
" \n",
" | \n",
" LLF | \n",
" AIC | \n",
" BIC | \n",
" chi2 | \n",
" p-value(chi2) | \n",
" n | \n",
" k | \n",
" dof | \n",
"
\n",
" \n",
" | scale | \n",
" | \n",
" | \n",
" | \n",
" | \n",
" | \n",
" | \n",
" | \n",
" | \n",
"
\n",
" \n",
" \n",
" \n",
" | plain | \n",
" -7.394196 | \n",
" 18.788392 | \n",
" 22.909278 | \n",
" 23.105937 | \n",
" 0.999971 | \n",
" 58 | \n",
" 2 | \n",
" 56 | \n",
"
\n",
" \n",
" | log | \n",
" -6.795379 | \n",
" 17.590757 | \n",
" 21.711643 | \n",
" 19.340467 | \n",
" 0.999999 | \n",
" 58 | \n",
" 2 | \n",
" 56 | \n",
"
\n",
" \n",
"
\n",
"