LinReg - 1var, 2param.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "03f52b3a",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import scipy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "88bd4c38",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "43ab148d",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "7e2e3bf7",
   "metadata": {},
   "outputs": [],
   "source": [
    "def SOS(A):\n",
    "    x=0\n",
    "    for i in A:\n",
    "        x+=i**2\n",
    "    return x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f3f48736",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "2061a88d",
   "metadata": {},
   "outputs": [],
   "source": [
    "#y=a+b*(x-x_bar)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "bae8644a",
   "metadata": {},
   "outputs": [],
   "source": [
    "X=np.array([2,2.2,2.4,2.4,2.6,2.6,2.8,3])\n",
    "Y=np.array([84,76,73,51,68,48,38,23])\n",
    "n=len(X)\n",
    "p=2                                       #edit manually"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "453f6bb7",
   "metadata": {},
   "outputs": [],
   "source": [
    "x_bar=np.mean(X)\n",
    "y_bar=np.mean(Y)\n",
    "a=y_bar\n",
    "x=0\n",
    "b=sum(((X-x_bar)*(Y-y_bar)))/SOS(X-x_bar)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "1277bb07",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "y_hat= 57.625 + (-59.30555555555555)*(x-2.5)\n"
     ]
    }
   ],
   "source": [
    "print(\"y_hat= {} + ({})*(x-{})\".format(a,b,x_bar))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "6cdf1b15",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([87.27777778, 75.41666667, 63.55555556, 63.55555556, 51.69444444,\n",
       "       51.69444444, 39.83333333, 27.97222222])"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "y_hat=a+b*(X-x_bar)\n",
    "y_hat"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "9f981229",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "565.5277777777779"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "SR=SOS(Y-y_hat)\n",
    "SR"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "85194bfe",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "94.25462962962966"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "s2=SR/(n-p)\n",
    "s2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "bffd9c48",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "9.70848235460258"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "s=np.sqrt(s2)\n",
    "s       #lack of fit error"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "b07dc554",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "s**2/sigma**2 = 377.01851851851865\n",
      "p-value = 0.0\n",
      "Reject H0\n"
     ]
    }
   ],
   "source": [
    "sigma=0.5                                       #given experimental error value, edit manually\n",
    "print(\"s**2/sigma**2 = {}\".format(s**2/sigma**2))\n",
    "p_val=1 - scipy.stats.chi2.cdf(s**2/sigma**2, (n-p))\n",
    "print(\"p-value = {}\".format(p_val))\n",
    "alpha=0.1                                      #given/chosen, edit manually\n",
    "if alpha==None:\n",
    "    alpha=1-p_val\n",
    "if p_val>=1-alpha:\n",
    "    print(\"H0 holds\")\n",
    "else:\n",
    "    print(\"Reject H0\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "14e1b09c",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "y_hat = -59.30555555555555*x+(205.88888888888889)\n"
     ]
    }
   ],
   "source": [
    "intercept=a-b*x_bar\n",
    "slope=b\n",
    "print(\"y_hat = {}*x+({})\".format(slope,intercept))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "1700b794",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Standard Error of intercept = 28.809102338879235\n",
      "Standard Error of co-efficient or slope(b) = 11.441556179949043\n"
     ]
    }
   ],
   "source": [
    "SE_intercept=s*((1/n)+x_bar**2/SOS(X-x_bar))**(0.5)\n",
    "SE_b=s/(np.sqrt(SOS(X-x_bar)))\n",
    "print(\"Standard Error of intercept = {}\".format(SE_intercept))\n",
    "print(\"Standard Error of co-efficient or slope(b) = {}\".format(SE_b))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "ccfd8539",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "t-ratio of intercept = 7.146661026332367\n",
      "t-ratio of co-efficient or slope(b) = -5.183346969836728\n"
     ]
    }
   ],
   "source": [
    "intercept_given=0                                #edit manually if provided\n",
    "b_given=0                                        #edit manually if provided\n",
    "t_ratio_intercept=(intercept-intercept_given)/SE_intercept\n",
    "t_ratio_b=(b-b_given)/SE_b\n",
    "print(\"t-ratio of intercept = {}\".format(t_ratio_intercept))\n",
    "print(\"t-ratio of co-efficient or slope(b) = {}\".format(t_ratio_b))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "ddbbbc16",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "p-value for intercept = 0.000189170391448501\n",
      "p-value for co-efficient or slope(b) = 0.9989761661075788\n"
     ]
    }
   ],
   "source": [
    "p_intercept=1-scipy.stats.t.cdf(t_ratio_intercept, (n-p))\n",
    "p_b=1-scipy.stats.t.cdf(t_ratio_b, (n-p))\n",
    "print(\"p-value for intercept = {}\".format(p_intercept))\n",
    "print(\"p-value for co-efficient or slope(b) = {}\".format(p_b))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "01fe5cfb",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Model: DOF= 2 SM= 29097.472222222223 SM2= 14548.736111111111\n",
      "Error: DOF= 6 SR= 565.5277777777779 SR2= 94.25462962962966\n",
      "Total: DOF= 8 ST= 29663.0\n",
      "F-ratio: 154.3556657989095\n"
     ]
    }
   ],
   "source": [
    "SM=SOS(y_hat)\n",
    "print(\"Model: DOF= {} SM= {} SM2= {}\".format(p,SM,SM/p))\n",
    "print(\"Error: DOF= {} SR= {} SR2= {}\".format(n-p,SR,SR/(n-p)))\n",
    "print(\"Total: DOF= {} ST= {}\".format(n,SM+SR))\n",
    "print(\"F-ratio: {}\".format((SM/p)/(SR/(n-p))))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "87e46998",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "p-value of F-ratio test: = 6.929734719940406e-06\n"
     ]
    }
   ],
   "source": [
    "p_F=1-scipy.stats.f.cdf((SM/p)/(SR/(n-p)),dfn=p,dfd=n-p)\n",
    "print(\"p-value of F-ratio test: = {}\".format(p_F))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1526e8cc",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c1a49bb8",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "468fd527",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a60000d6",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "39d3e5d5",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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",
   "version": "3.9.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}