{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Statystyka i Analiza danych\n",
    "# Laboratorium 8 - Korelacja i regresja (część 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Ćwiczenie 1: Wyznaczanie korelacji\n",
    "\n",
    "W ćwiczeniu tym sprawdzimy korelację między dwoma sygnałami. W szczególności sprawdzimy, jak zmiany w definicji oraz we wzajemnym położeniu sygnałów wpływają na wartość korelacji.\n",
    "\n",
    "Zaczniemy od przygotowania danych. Wykonaj poniższy kod, aby przygotować dane do ćwiczenia."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "t <- seq(0,4,0.2)\n",
    "n <- length(t)\n",
    "\n",
    "prepareData <- function(x1, x2){\n",
    "  return(data.frame(t, x1, x2))\n",
    "}\n",
    "# Tworzymy różne pary sygnałów\n",
    "sygnaly_1 <- prepareData(x <- runif(n),1.5*x+runif(n)/2)\n",
    "sygnaly_2 <- prepareData(sin(t), cos(t))\n",
    "sygnaly_3 <- prepareData(runif(n), runif(n))\n",
    "sygnaly_4 <- prepareData(runif(n), 10*runif(n)-5)\n",
    "sygnaly_5 <- prepareData(x <- runif(n), c(x[n],x[1:n-1]))\n",
    "sygnaly_6 <- prepareData(sin(t), 5*sin(t))\n",
    "sygnaly_7 <- prepareData(5*sin(t)+rnorm(n), 5*sin(t)+rnorm(n))\n",
    "\n",
    "# Usuwamy niepotrzebne zmienne (w tym funkcje)\n",
    "rm(n, t, x, prepareData)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Każda ze stworzonych wcześniej ramek danych zawiera informacje o wartościach dwóch sygnałów, `x1` i `x2` w czasie `t`. Wyświetl zawartość ramki `sygnaly_1`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Przypisz wektory `sygnaly_1$t`, `sygnaly_1$x1` i `sygnaly_1$x2` do zmiennych `t`, `x1` i `x2`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Korzystając z funkcji `plot` wykonaj wykres rozrzutu (*scattered plot*), nazwij osie *X1* i *X2* i nadaj wykresowi tytuł *Sygnały 1*."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plot(x2 ~ x1, ...) - reszta parametrów do uzupełnienia"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Powstaje pytanie czy wartości sygnałów są ze sobą powiązane - innymi słowy czy są skorelowane? Aby to sprawdzić możemy policzyć pewien punkt odniesienia (wartość średnią) i dla każdego `t` sprawdzić, po której stronie wartości średniej znajduje się wartość sygnału. Jeżeli dla dwóch sygnałów znak będzie taki sam, to ich iloczyn będzie dodatni, w przeciwnym przypadku ujemny."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Wykonaj poniższy kod, aby zdefiniować funkcje `signalPlotWithMeans` oraz `scatterPlotWithMeans` pozwalające na umieszczanie na wykresie informacji o przebiegach i wartościach średnich sygnałów."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "signalPlotWithMeans <- function(t, x1, x2, title){\n",
    "  min <- min(x1, x2)\n",
    "  max <- max(x1,x2)\n",
    "  plot(x1 ~ t, ylim=c(min, max), type=\"o\", col=\"red\", ylab=\"value\", main=title)\n",
    "  abline(h=mean(x1), col=\"red\")\n",
    "  lines(x2 ~ t, col=\"green\",type=\"o\")\n",
    "  abline(h=mean(x2), col=\"green\")\n",
    "  legend('topright', c(\"x1\", \"x2\", \"E(x1)\", \"E(x2)\"), lty=1, col=c('red', 'green', 'red',' green'), bty='n', cex=.75)\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "scatterPlotWithMeans <- function(x1, x2, title){\n",
    "    plot(x2~x1, xlab=\"x1\", ylab=\"x2\", main=title)\n",
    "    abline(v = mean(x1), col=\"red\", lty=3)\n",
    "    abline(h = mean(x2), col=\"green\", lty=3)\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Korzystając z funkcji `plotWithMeans` stwórz wykres wartości sygnału w zależności od czasu z naniesionymi wartościami średnimi sygnałów."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "signalPlotWithMeans(t, x1, x2, \"sygnaly_1\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "scatterPlotWithMeans(x1, x2, \"sygnaly_1\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Oblicz wartości średnie sygnałów `x1` i `x2`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Oblicz iloczyny odchyleń wartości sygnałów od średnich `(x1-ex1)*(x2-ex2)`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Uśredniając takie iloczyny otrzymujemy kowariancję dwóch sygnałów. Oblicz kowariancję (z próby) dla `x1` i `x2`. Zweryfikuj poprawność obliczeń korzystając z wbudowanej funkcji cov."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Kowariancja jest zależna od wartości bezwględnych. Jej podzielenie przez iloczyn odchyleń standardowych normuje wartość współczynnika do przedziału od -1 do 1 - jest to korelacja. Oblicz wartość współczynnika korelacji $\\frac{cov(x1, x2)}{sd(x1)\\cdot sd(x2)}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Na podstawie wykonanych wcześniej poleceń uzupełnij funkcję `analyzeCorr`, która przyjmie jako argumenty `t`, `x1`, `x2` i tytuł (`title`), stworzy wykres rozrzutu i wykres zależności sygnałów od czasu wraz z średnimi (wykresy powinny być wyświetlone obok siebie) i zwróci wartość współczynnika korelacji dla `x1` i `x2`. Do wykreślenia wykresów obok siebie możesz skorzystać z polecenia `par(mfrow=c(1,2))`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "analyzeCor <- function(t, x1, x2, title){\n",
    "    # Ustawiamy dwa wykresy obok siebie\n",
    "    prev  <- par(mfrow=c(1,2))\n",
    "     # 1. wykres scatterPlotWithMeans \n",
    "    scatterPlotWithMeans(x1, x2, title)\n",
    "    # 2. wkres signalPlotWithMeans\n",
    "    signalPlotWithMeans(t, x1, x2, title)\n",
    "    # Przywracamy poprzednie ustawienia\n",
    "    par(prev)\n",
    "    # Wyznaczamy kowariancję i korelację, zwracamy tę ostatnią\n",
    "    # ... do uzupełnienia\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Sprawdzić działanie funkcji `analyzeCorr`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cor  <- analyzeCor(t, x1, x2, \"sygnały 1\")\n",
    "cor"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Za pomocą funkcji `analyzeCorr` przeanalizuj wszystkie pary sygnałów wykorzystując poniższą pętlę."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cor  <- vector(length = 7)\n",
    "names(cor)  <- paste(\"sygnaly\", 1:7, sep=\"_\")\n",
    "for (name in names(cor)){\n",
    "  signals <- eval(as.name(name))\n",
    "  cor[name] <- # ... do uzupełnienia\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A następnie wyświetl współczynniki korelacji dla poszczególnych par sygnałów (zmnienna `cor`)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cor"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Porównaj powyższe wyniki z wartościami uzykanymi za pomocą funkcji wbudowanej `cor`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cor  <- vector(length = 7)\n",
    "names(cor)  <- paste(\"sygnaly\", 1:7, sep=\"_\")\n",
    "for (name in names(cor)){\n",
    "  signals <- eval(as.name(name))\n",
    "  cor[name] <- # ... do uzupełnienia\n",
    "}\n",
    "cor"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Ćwiczenie 2: Wyznaczanie współczynników modelu regresji\n",
    "\n",
    "W ćwiczeniu tym wyznaczymy współczynniki w modelu regresji stosując metodę najmnieszych kwadratów.\n",
    "\n",
    "Podobnie jak poprzednio, zaczniemy od przygotowania danych. Wykonaj poniższy kod, aby przygotować dane do ćwiczenia."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "regresja <- read.csv(url(\"http://www.cs.put.poznan.pl/swilk/siad/8-cw2.csv\"), sep=\";\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ramka zawiera trzy kolumny, `x`, `yp=f(x)`  i `y=f(x)+epsilon`. Postaramy się wyestymować parametry `a` i `b`, tak aby jak najlepiej przybliżyć wartość `y`. Pierwsze podejście będzie polegać na ręcznym dopasowaniu wartości `a^` i `b^`, tak, aby suma różnic kwadratów była jak najmniejsza. \n",
    "\n",
    "Uzupełnij definicję funkcji `plotReg`, która wykreśla `yp`, `y` oraz `y^=a^*x+b^` i zwraca wartość sumy kwadratów różnic `(y^-y)^2`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "plotReg <- function(reg, at, bt){\n",
    "    # Oblicz y^=a^*x+b^ (yt)\n",
    "    # ... do uzupełnenia\n",
    "    # Oblicz sumę (y^-y)^2:\n",
    "    # ... do uzupełnienia\n",
    "  \n",
    "    plot(reg$yp ~ reg$x, t=\"l\", col=\"red\", xlab=\"x\", ylab=\"\")\n",
    "    points(reg$y ~ reg$x, col=\"blue\")\n",
    "    lines(yt ~ reg$x, t=\"l\", col=\"green\")\n",
    "  \n",
    "    # Funkcja powinna zwracać sumę kwadratów różnic:\n",
    "    # ... do uzupełnienia\n",
    "}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Spróbuj dopasować wartości `a^` i `b^` korzystając ze zdefiniowanej powyżej funkcji."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Następnie wykorzystamy funkcję `lm`. Pierwszym argumentem tej metody jest `formula`, którą zdefinijemy jako `y~x`, drugim jest `data=regresja`. Otrzymany model regresji zapisz do zmiennej `model`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Korzystając ze zmiennej model odczytaj 'a^' i 'b^'."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Obszerniejsze wyniki można uzyskać korzystając z funkcji `summary`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Ćwiczenie 3: Testowanie istotności statystycznej korelacji\n",
    "\n",
    "W tym ćwiczeniu sprawdzimy, czy korelacja między zmiennymi jest statystycznie istotna. W tym celu skorzystamy ze statystyki *t*.\n",
    "\n",
    "Zwyczajowo już zaczynamy od załadowania danych wykonyjąc poniższy kod."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dzieci <- read.csv(url(\"http://www.cs.put.poznan.pl/swilk/siad/8-cw3.csv\"), sep=\";\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Dane te zawierają informację o wieku w latach (`x`) i wzroście w cm (`y`) dla 15-osobowej grupy wylosowanej z populacji dzieci i młodzieży. Wyświetl zawartości ramki `dzieci`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Korzystajac z funkcji `plot` utwórz wykres rozrzutu dla tych danych. Opisz odpowiednio osie wykresu."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Oblicz średnie i odchylenia standardowe obu zmiennych oraz liczbę obserwacji."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Oblicz samodzielnie współczynnik korelacji korzystając ze wzoru $\\frac{cov(x, y)}{sd(x) \\cdot sd(y)}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Oblicz także współczynnik korelacji korzystając z funkcji `cor`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Aby zweryfikować hipotezę mówiącą, że współczynnik korelacji jest niezerowy (test dwustronny), oblicz wartość statystyki $t=\\frac{r}{\\sqrt{1-r^2}}\\cdot \\sqrt{n-2}$, gdzie $r$ to współczynnik korelacji dla próby."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Oblicz p-wartość dla wyznaczonej statystyki pamiętając, że mamy do czynienia z testem dwustronnym. Jaki będzie wynik testowania hipotez (zakładając poziom istotności równy 0.05)?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Wyznacz model regresji korzystając z funkcji `lm`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Utwórz wykres modelu regresji korzystając z funkcji `abline` (wcześniej narysuj ponownie wykres rozrzutu). W jakich warunkach uzyskany model może mieć zastosowanie?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "R",
   "language": "R",
   "name": "ir"
  },
  "language_info": {
   "codemirror_mode": "r",
   "file_extension": ".r",
   "mimetype": "text/x-r-source",
   "name": "R",
   "pygments_lexer": "r",
   "version": "3.6.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}