code/mathematica/wolfram.ipynb

{
  "cells": [
    {
      "cell_type": "markdown",
      "source": "[矢吹太朗『コンピュータでとく数学』（オーム社, 2024）](https://github.com/taroyabuki/comath)\n\n\n",
      "metadata": {}
    },
    {
      "cell_type": "markdown",
      "source": "# 1 実行環境",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "data = {{1, 7}, {3, 1}, {6, 6}, {10, 14}};\n",
        "model = LinearModelFit[data, X1, X1]\n",
        "model[\"BestFitParameters\"]"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 2 数と変数",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "2 * (-3)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "2 (-3)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "(1 + 2) * 3"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "2^10"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "-2 < -1"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "2 + 2 == 5"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "If[7 < 5, 10, 20]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x<1"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x == y"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x^2 - 1 == (x + 1) (x - 1) // Simplify"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Not[0 < 1] (* 方法1 *)\n",
        "! (0 < 1)  (* 方法2 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Or[0 < 1, 2 > 3]   (* 方法1 *)\n",
        "(0 < 1) || (2 > 3) (* 方法2 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "And[0 < 1, 2 > 3]  (* 方法1 *)\n",
        "(0 < 1) && (2 > 3) (* 方法2 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Not[10<x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x = 5; x == 5"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "a = 1 + 2;\n",
        "b = 9;\n",
        "a(b + 1)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "a = 1 + 2; b = 9; a * (b + 1)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "a = 1 + 2"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "a = 3;\n",
        "Clear[a]; (* 変数を記号にする．*)\n",
        "Expand[(a + 1)^2]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x1 = 2; x2 = 3; x1 + x2"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Subscript[x, 1] = 2; Subscript[x, 2] = 3; Subscript[x, 1] + Subscript[x, 2]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x = 1; y = x + 1; x = 2; y"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x = 1;\n",
        "y := x + 1; (* yは「2」ではなく「x + 1」になる．*)\n",
        "x = 2;\n",
        "y           (* 「x + 1」は「2 + 1」つまり3．*)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f = 2 x + 3;\n",
        "f /. x -> 5"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "g = a + b;\n",
        "g /. {a -> x, b -> y}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f = Function[x, 2 x + 3];\n",
        "f[5]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[f];\n",
        "f[x_] := 2 x + 3\n",
        "f[5]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[f, a];\n",
        "f = Function[x, 2 x + 3];\n",
        "g = f[a];\n",
        "{f[5], g /. a -> 5}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[f];\n",
        "f[x_] := 1 / x\n",
        "f[1]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f1[x_] := Piecewise[{{1/x, x != 0}}, Undefined]\n",
        "\n",
        "f2[0] = Undefined;\n",
        "f2[x_] := 1/x\n",
        "\n",
        "f3[0] = Undefined;\n",
        "f3[x_ /; x != 0] := 1/x\n",
        "\n",
        "f4[x_] := If[x != 0, 1/x, Undefined]\n",
        "\n",
        "f5[x_] := Which[x != 0, 1/x, True, Undefined]\n",
        "\n",
        "{f1[1], f2[1], f3[1], f4[1], f5[1]} (* 全て1 *)\n",
        "{f1[0], f2[0], f3[0], f4[0], f5[0]} (* 全てUndefined *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Function[x, 2 x + 3][5]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[f];\n",
        "f[x_, y_] := x + y\n",
        "f[2, 3]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[g];\n",
        "g[x_] := x[[1]] + x[[2]]\n",
        "x = {2, 3}; g[x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "g[{x1_, x2_}] := x1 + x2\n",
        "g[x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Apply[f, x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "g[{2, 3}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Expand[(x+1)^2]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "N[Sqrt[2],30]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "pi2 = FromDigits[RealDigits[N[Pi], 2], 2]\n",
        "pi10 = FromDigits[RealDigits[N[Pi], 10], 10]\n",
        "Abs[Pi - pi2] < Abs[Pi - pi10] (* True *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "0.1 + 0.2 == 0.3"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Chop[0.1 + 0.2 - 0.3] == 0"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "1/10+2/10==3/10"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "(*「‘」はシングルクォートではなくバッククォート *)\n",
        "Block[{Internal`$EqualTolerance = 0.}, 0.1 + 0.2 == 0.3] (* False *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Chop[0.1 + 0.2 - 0.3] == 0 (* True *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[x];\n",
        "Simplify[Sin[x]^2 + Cos[x]^2]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "FullSimplify[Sqrt[5 + 2 Sqrt[6]]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Simplify[Sqrt[(x - 1)^2], x - 1 >= 0]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[Sqrt[(x - 1)^2] == x - 1, x, Reals]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[a, b];\n",
        "Reduce[Sqrt[a] Sqrt[b] == Sqrt[a b], Reals]"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 3 データ構造",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "v = {2, 3, 5}; Length[v]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "v[[3]] = 0.5; v"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Range[5]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Range[0, 1, 0.1]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Subdivide[0, 100, 4]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "v = {2, 3};\n",
        "1.1 v"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "u = {10, 20}; v = {2, 3};\n",
        "u + v"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "v + 1"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "u = {10, 20}; v = {2, 3};\n",
        "u . v"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "a = {2, 3, 4}; b = a; b[[3]] = 0.5; a"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "v = {2, -1, 3, -2};\n",
        "Cases[v, x_ /; x > 0]         (* パターンマッチングによる抽出 *)\n",
        "Select[v, Function[x, x > 0]] (* 関数による抽出 *)\n",
        "Select[v, Positive]           (* 組込み関数の利用 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "v = {2, -1, 3, -2};\n",
        "UnitStep[v]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "v = {2, -1, 3, -2};\n",
        "n = Length[v];      (* vのサイズ *)\n",
        "u = Table[Null, n]; (* Nullは「値がない」ということ．*)\n",
        "Do[u[[i]] = If[v[[i]] < 0, 0, 1], {i, 1, n}];\n",
        "u"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Table[If[x < 0, 0, 1], {x, v}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "v = {2, -1, 3, -2};\n",
        "f = Function[x, If[x < 0, 0, 1]];\n",
        "Map[f, v]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "v = {2, -1, 3, -2};\n",
        "f = Function[x, If[x < 0, 0, 1], Listable];\n",
        "f[v]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "u = {1, 7, 2, 9}; v = {2, 3, 5, 7};\n",
        "f = Function[{a, b}, If[a < b, -1, 1]];\n",
        "MapThread[f, {u, v}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x = <|\"apple\" -> \"りんご\", \"orange\" -> \"みかん\"|>;\n",
        "x[\"orange\"]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "AppendTo[x, \"grape\" -> \"ぶどう\"];\n",
        "x[\"grape\"]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x[\"apple\"] =.\n",
        "KeyExistsQ[x, \"apple\"]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[x];\n",
        "x[\"apple\"] = \"りんご\";\n",
        "x[\"orange\"] = \"みかん\";\n",
        "\n",
        "x[\"orange\"]            (* みかん *)\n",
        "\n",
        "x[\"grape\"] = \"ぶどう\";\n",
        "x[\"grape\"]             (* ぶどう *)\n",
        "\n",
        "x[\"apple\"] =.\n",
        "Head[x[\"apple\"]] =!= x (* False *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "df = Transpose[Dataset[<|\"name\" -> {\"A\", \"B\", \"C\"},\n",
        "                         \"english\" -> {60, 90, 70},\n",
        "                         \"math\" -> {70, 80, 90},\n",
        "                         \"gender\" -> {\"f\", \"m\", \"m\"}|>]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "df = Dataset[{\n",
        "   <|\"name\" -> \"A\", \"english\" -> 60, \"math\" -> 70, \"gender\" -> \"f\"|>,\n",
        "   <|\"name\" -> \"B\", \"english\" -> 90, \"math\" -> 80, \"gender\" -> \"m\"|>,\n",
        "   <|\"name\" -> \"C\", \"english\" -> 70, \"math\" -> 90, \"gender\" -> \"m\"|>}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "df[All, {\"english\", \"math\"}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Normal[df[[All, \"english\"]]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "m = Values[Normal[df[[All, {2, 3}]]]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "{english, math} = Transpose[m]"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 4 可視化と方程式",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Plot[x^2 + 2 x - 4, {x, -5, 3}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x = Subdivide[-5, 3, 100];\n",
        "y = x^2 + 2 x - 4;\n",
        "ListLinePlot[Transpose[{x, y}]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Plot[Sin[102 x], {x, 0, 2 Pi}, PlotPoints -> 100]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Plot3D[x^2 + y^2, {x, -1, 1}, {y, -1, 1}, AxesLabel -> {\"x\", \"y\"}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "ContourPlot[x^2 + y^2, {x, -1, 1}, {y, -1, 1}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[x, y];\n",
        "reg1 = ImplicitRegion[x^2 + y^2 == 1, {x, y}];\n",
        "Region[reg1, Axes -> True]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "reg2 = ImplicitRegion[x^2 + y^2 <= 1, {x, y}];\n",
        "Region[reg2, Axes -> True]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "reg = ImplicitRegion[Or[2 x + 3 y == 8, 5 x - 7 y == -9], {x, y}];\n",
        "Region[reg, Axes -> True]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "reg = ImplicitRegion[And[y <= x, y >= x^2], {x, y}];\n",
        "RegionPlot[reg]                                                    (* ① *)\n",
        "Plot[{x, x^2}, {x, 0, 1}, Filling -> {1 -> {2}}, AspectRatio -> 1] (* ② *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "{RegionMeasure[reg1], RegionMeasure[reg2]}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "RegionMeasure[reg]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[x];\n",
        "{a, b} = Sort[SolveValues[{x == x^2}, x]];\n",
        "Integrate[x - x^2, {x, a, b}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Manipulate[\n",
        " ContourPlot[x^3 + y^3 - 3 x y == c, {x, -2, 2}, {y, -2, 2}],\n",
        " {{c, 0}, -1, 1}] (* cは-1以上1以下で，初期値は0 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "SolveValues[x^2 + 2 x - 4 == 0, x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "{a, b} = SolveValues[x^2 + 2 x - 4 == 0, x]\n",
        "a + b\n",
        "\n",
        "tmp = Solve[x^2 + 2 x - 4 == 0, x]\n",
        "{a, b} = x /. tmp;\n",
        "a + b\n",
        "\n",
        "tmp = Reduce[x^2 + 2 x - 4 == 0, x]\n",
        "{a, b} = x /. {ToRules[tmp]};\n",
        "a + b"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "n = 3; Simplify[Total[SolveValues[x^n + 2 x - 4 == 0, x]]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "sol = SolveValues[{2 x + 3 y == 8, 5 x - 7 y == -9}, {x, y}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "{{x1, y1}} = sol; x1 + y1"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[x_] := 2^x + Sin[x]\n",
        "FindRoot[f[x] == 0, {x, 0}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[x^2 + 2 x - 4 < 0, x]"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 5 論理式",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "expr = Exists[x, Element[x, Reals], x^2 == 2];\n",
        "Reduce[expr]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[Implies[x > 10, x > 11]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[ForAll[x, Element[x, Reals], Implies[x > 10, x > 11]]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "BooleanConvert[Implies[A, B], \"OR\"] (* 含意 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "BooleanConvert[And[A, B], \"OR\"] (* 論理積 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "{BooleanConvert[Not[A]] == BooleanConvert[Nand[A, A]],\n",
        " BooleanConvert[Or[A, B]] == BooleanConvert[Nand[Not[A], Not[B]]]}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[Exists[x, Element[x, Reals], x^2 == 2]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[Exists[x, x^2 == 2], Reals]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[Exists[x, Element[x, Rationals], x^2 == 2]] (* False *)\n",
        "Reduce[Exists[x, x^2 == 2], Rationals]             (* False *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "expr = ForAll[b, Element[b, Reals], Exists[n, Element[n, Integers], n > b]];\n",
        "Reduce[expr]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "expr1 = ForAll[b,\n",
        "   Element[b, Reals], Exists[n, And[Element[n, Integers], n > b]]];\n",
        "Reduce[expr1] (* True *)\n",
        "\n",
        "expr2 = ForAll[b,\n",
        "   Implies[Element[b, Reals], Exists[n, Element[n, Integers], n > b]]];\n",
        "Reduce[expr2] (* 失敗 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[Exists[x, a x + b == 0]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[Exists[x, Element[x, Reals], x^2 + a^2 == 0]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[Exists[x, Element[x, Reals], x^2 + a^2 < 0]] (* False *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[Exists[x, Element[x, Reals], x^2 + a^2 < 0], Complexes] // Simplify"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[Not[Exists[{n, a, b, c}, And[n >= 3, a^n + b^n == c^n]]],\n",
        " PositiveIntegers]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[Not[Exists[{a, b, c}, a^4 + b^4 == c^4]], PositiveIntegers]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[Not[Exists[{a, b, c}, a^4 + b^4 == c^2]], PositiveIntegers] (* 失敗 *)"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 6 1次元のデータ",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "a = {36, 43, 53, 55, 56, 56, 57, 60, 61, 73};\n",
        "b = {34, 39, 39, 49, 50, 52, 52, 55, 83, 97};\n",
        "Histogram[a]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "HistogramList[a, {20, 80, 20}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x = {7, 3, 1, 3, 4, 7, 7, 7, 10, 3};\n",
        "f = Counts[x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "BoxWhiskerChart[{a, b}, \"Outliers\", ChartLabels -> {\"A\", \"B\"}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "a = {36, 43, 53, 55, 56, 56, 57, 60, 61, 73};\n",
        "Mean[a]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "b = {34, 39, 39, 49, 50, 52, 52, 55, 83, 97};\n",
        "Total[b]/Length[b]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Mean[a - Mean[a]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Variance[a]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Total[(b - Mean[b])^2]/(Length[b] - 1) // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "z = Standardize[a]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "{Mean[z], StandardDeviation[z]}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "(a - Mean[a])/StandardDeviation[a]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "StandardDeviation[a] z + Mean[a]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "10 * z + 50"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 7 2次元のデータ",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x = {35, 45, 55, 65, 75}; y = {114, 124, 143, 158, 166};\n",
        "ListPlot[Transpose[{x, y}]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x = {35, 45, 55, 65, 75}; y = {114, 124, 143, 158, 166};\n",
        "Covariance[x, y]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Covariance[Transpose[{x, y}]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "(x - Mean[x]) . (y - Mean[y])/(Length[x] - 1)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Correlation[x, y] // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x = {35, 45, 55, 65, 75}; y = {114, 124, 143, 158, 166};\n",
        "data = Thread[{x, y}]; (* x, yを列とする行列 *)\n",
        "model = LinearModelFit[data, X, X]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "model[40]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Show[ListPlot[data], Plot[model[x], {x, 35, 75}]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "L = Total[(y - (a x + b))^2]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "sol = SolveAlways[L == p (a - q)^2 + r (b - (s a + t))^2 + u, {a, b}]\n",
        "{q, s q + t} /. sol[[1]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "a = Covariance[x, y]/Variance[x]; b = Mean[y] - a Mean[x];\n",
        "{a, b} // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "anscombe = ExampleData[{\"Statistics\", \"AnscombeRegressionLines\"}];\n",
        "x1 = anscombe[[All, 1]]; y1 = anscombe[[All, 5]]; data = Thread[{x1, y1}];\n",
        "Correlation[x1, y1]\n",
        "model = LinearModelFit[data, X, X]\n",
        "Show[ListPlot[data], Plot[model[x], {x, 0, 21}]]"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 8 確率変数と確率分布",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = DiscreteUniformDistribution[{1, 6}];\n",
        "PDF[dist][2]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Probability[X == 2, Distributed[X, dist]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "data = RandomVariate[dist, 1000];\n",
        "Histogram[data] (* 結果は割愛 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Show[Histogram[data, {0.5, 6.5, 1}, \"PDF\"],\n",
        " DiscretePlot[PDF[dist][x], {x, 1, 6}]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = BernoulliDistribution[3/10];\n",
        "data = RandomVariate[dist, 1000];\n",
        "Counts[data]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = BinomialDistribution[10, 3/10];\n",
        "PDF[dist][3]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Probability[X == 3, Distributed[X, dist]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = BinomialDistribution[n, p];\n",
        "PDF[dist]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "n = 10; p = 3/10; dist = BinomialDistribution[n, p];\n",
        "data = RandomVariate[dist, 1000];\n",
        "Histogram[data] (* 結果は割愛 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Show[Histogram[data, {-0.5, n + 0.5, 1}, \"PDF\"],\n",
        " DiscretePlot[PDF[dist][x], {x, 0, n}]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = BinomialDistribution[10, 3/10];\n",
        "CDF[dist][3]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Probability[X <= 3, Distributed[X, dist]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Sum[PDF[dist][k], {k, 0, 3}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = UniformDistribution[{0, 360}];\n",
        "{CDF[dist][200], CDF[dist][150], CDF[dist][200] - CDF[dist][150]}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Probability[150 <= X <= 200, Distributed[X, dist]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Integrate[PDF[dist][x], {x, 150, 200}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Integrate[PDF[dist][t], {t, 0, x},\n",
        " Assumptions -> Element[x, Reals]] (* xは実数と仮定する．*)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "D[x/360, x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "data = RandomVariate[dist, 1000];\n",
        "Histogram[data] (* 結果は割愛 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "data = RandomVariate[dist, 1000];\n",
        "Show[Histogram[data, {\"Raw\", \"Sturges\"}, \"PDF\"],\n",
        " Plot[PDF[dist][x], {x, 0, 360}]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = NormalDistribution[6, 2];\n",
        "CDF[dist][6 + 3 2] - CDF[dist][6 - 3 2] // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Probability[6 - 3 2 <= X <= 6 + 3 2, Distributed[X, dist]] // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Integrate[PDF[dist][x], {x, 6 - 3 2, 6 + 3 2}] // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[mu, sigma, x];\n",
        "dist = NormalDistribution[mu, sigma];\n",
        "{a, b} = {mu - 3 sigma, mu + 3 sigma};\n",
        "CDF[dist][b] - CDF[dist][a] // N                    (* 方法1 *)\n",
        "Probability[a <= X <= b, Distributed[X, dist]] // N (* 方法2 *)\n",
        "Integrate[PDF[dist][x], {x, a, b}] // N             (* 方法3 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = NormalDistribution[mu, sigma];\n",
        "PDF[dist][x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Plot[{PDF[NormalDistribution[0, 1]][x],\n",
        "      PDF[NormalDistribution[2, 1]][x]}, {x, -6, 6}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Plot[{PDF[NormalDistribution[0, 1]][x],\n",
        "      PDF[NormalDistribution[0, 2]][x]}, {x, -6, 6}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Xs = {0, 100, 1000, 10000}; Ps = {0.9, 0.08, 0.015, 0.005};\n",
        "tmp = Piecewise[Thread[{Ps, Thread[x == Xs]}]];\n",
        "dist = ProbabilityDistribution[tmp, {x, 0, 10000, 1}, (* 確率分布の定義 *)\n",
        "   Method -> \"Normalize\"]; (* 念のため合計を1にする．*)\n",
        "data = RandomVariate[dist, 1000];\n",
        "Counts[data]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = ProbabilityDistribution[Abs[x], {x, -2, 2}, (* 確率分布の定義 *)\n",
        "   Method -> \"Normalize\"]; (*“全確率”が1にならない場合の備え *)\n",
        "data = RandomVariate[dist, 1000];\n",
        "Histogram[data]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[f, F, t, x];\n",
        "f[x_] := RealAbs[x]                            (* 手順1 *)\n",
        "F[x_] := Evaluate[Integrate[f[t], {t, -1, x}]] (* 手順2 *)\n",
        "Finv = InverseFunction[F];                     (* 手順3 *)\n",
        "data = Table[Finv[RandomReal[]], {1000}];      (* 手順4 *)\n",
        "Histogram[data]                                (* 手順5 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "distY = UniformDistribution[{0, 1}];\n",
        "distX = TransformedDistribution[\n",
        "   Piecewise[{{-Sqrt[1 - 2 Y], Y <= 1/2}}, Sqrt[-1 + 2 Y]],\n",
        "   Distributed[Y, distY]];\n",
        "PDF[distX]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "distX = DiscreteUniformDistribution[{1, 6}];\n",
        "distY = TransformedDistribution[Mod[X^3, 4], Distributed[X, distX]];\n",
        "data = RandomVariate[distY, 1000];\n",
        "Histogram[data]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "distX = UniformDistribution[{0, 1}];\n",
        "distY = TransformedDistribution[X^2, Distributed[X, distX]];\n",
        "data = RandomVariate[distY, 1000];\n",
        "Histogram[data]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "PDF[distY]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "distX = ProbabilityDistribution[Abs[x], {x, -1, 1}];\n",
        "distY = TransformedDistribution[X^2, Distributed[X, distX]];\n",
        "PDF[distY]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = NormalDistribution[mu, sigma]; Clear[a, b];\n",
        "TransformedDistribution[a X + b, Distributed[X, dist]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Xs = {0, 100, 1000, 10000}; Ps = {0.9, 0.08, 0.015, 0.005};\n",
        "tmp = Piecewise[Thread[{Ps, Thread[x == Xs]}]];\n",
        "dist = ProbabilityDistribution[tmp, {x, 0, 10000, 1}];\n",
        "Expectation[X, Distributed[X, dist]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Mean[dist]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Sum[x PDF[dist][x], {x, Xs}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Xs . Ps"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Mean[RandomVariate[dist, 500000]] // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[n, p];\n",
        "dist = BinomialDistribution[n, p];\n",
        "Expectation[X, Distributed[X, dist]]       (* 方法1 *)\n",
        "Mean[dist]                                 (* 方法2 *)\n",
        "Sum[x PDF[dist][x], {x, 0, n}] // Simplify (* 方法3 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = ProbabilityDistribution[Abs[x], {x, -1, 1}];\n",
        "Integrate[x PDF[dist][x], {x, -1, 1}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Xs = {0, 100, 1000, 10000}; Ps = {0.9, 0.08, 0.015, 0.005};\n",
        "tmp = Piecewise[Thread[{Ps, Thread[x == Xs]}]];\n",
        "dist = ProbabilityDistribution[tmp, {x, 0, 10000, 1}];\n",
        "Variance[dist]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Expectation[(X - Mean[dist])^2, Distributed[X, dist]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Sum[(x - Mean[dist])^2 PDF[dist][x], {x, Xs}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "((Xs - Xs . Ps)^2) . Ps"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[n, p];\n",
        "dist = BinomialDistribution[n, p];\n",
        "Variance[dist]                                              (* 方法1 *)\n",
        "Expectation[(X - Mean[dist])^2, Distributed[X, dist]]       (* 方法2 *)\n",
        "Sum[(x - Mean[dist])^2 PDF[dist][x], {x, 0, n}] // Simplify (* 方法3 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = ProbabilityDistribution[Abs[x], {x, -1, 1}];\n",
        "Integrate[(x - Mean[dist])^2 PDF[dist][x], {x, -1, 1}]"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 9 多次元の確率分布",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "distX = DiscreteUniformDistribution[{1, 6}];\n",
        "dist = TransformedDistribution[{Max[X1, X2], Min[X1, X2]},\n",
        "   {Distributed[X1, distX], Distributed[X2, distX]}];\n",
        "probs = Table[{\n",
        "   Probability[{X, Y} == {x, y}, Distributed[{X, Y}, dist]], (* 確率 *)\n",
        "   {X, Y} == {x, y}},                                        (* 条件 *)\n",
        "  {x, 1, 6}, {y, 1, 6}]\n",
        "\n",
        "dist = ProbabilityDistribution[Piecewise[Flatten[probs, 1]], (* 作り直し *)\n",
        "   {X, 1, 6, 1}, {Y, 1, 6, 1}];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "PDF[MarginalDistribution[dist, 1]][x] // Simplify\n",
        "PDF[MarginalDistribution[dist, 2]][y] // Simplify"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Table[CDF[dist][{x, y}], {y, 1, 6}, {x, 1, 6}] // TableForm"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "c = Counts[Flatten[Table[{Max[x, y], Min[x, y]}, {x, 1, 6}, {y, 1, 6}], 1]]/36;\n",
        "dist = ProbabilityDistribution[Piecewise[KeyValueMap[{#2, {X, Y} == #1} &, c]],\n",
        "   {X, 1, 6, 1}, {Y, 1, 6, 1}];\n",
        "\n",
        "Mean[dist]                (* 平均 *)\n",
        "Variance[dist]            (* 分散 *)\n",
        "StandardDeviation[dist]   (* 標準偏差 *)\n",
        "Covariance[dist][[1, 2]]  (* 共分散 *)\n",
        "Correlation[dist][[1, 2]] (* 相関係数 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "{uX, uY} = Mean[dist]; {sX, sY} = StandardDeviation[dist];\n",
        "Expectation[{X, Y,           (* 平均 *)\n",
        "   (X - uX)^2, (Y - uY)^2,   (* 分散 *)\n",
        "   (X - uX) (Y - uY),        (* 共分散 *)\n",
        "   (X - uX) (Y - uY)/sX/sY}, (* 相関係数 *)\n",
        "  Distributed[{X, Y}, dist]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Sum[x Probability[X == x, Distributed[{X, Y}, dist]], {x, 1, 6}] (* 平均 *)\n",
        "Sum[(x - uX) (y - uY) PDF[dist][{x, y}], {x, 1, 6}, {y, 1, 6}]   (* 共分散 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = DiscreteUniformDistribution[{1, 6}];\n",
        "Probability[Conditioned[X == 2, X <= 3], Distributed[X, dist]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "c = Counts[Flatten[Table[{Max[x, y], Min[x, y]}, {x, 1, 6}, {y, 1, 6}], 1]]/36;\n",
        "dist = ProbabilityDistribution[Piecewise[KeyValueMap[{#2, {X, Y} == #1} &, c]],\n",
        "   {X, 1, 6, 1}, {Y, 1, 6, 1}];\n",
        "rule = Distributed[{X, Y}, dist];\n",
        "\n",
        "Table[Probability[Conditioned[X == x, Y == 3], rule], {x, 1, 6}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Table[\n",
        " Probability[And[X == x, Y == 3], rule]/Probability[Y == 3, rule], {x, 1, 6}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Expectation[Conditioned[X, Y == 3], rule]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Sum[x Probability[Conditioned[X == x, Y == 3], rule], {x, 1, 6}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Table[Probability[And[X <= x, Y <= y], rule], {x, 1, 6}, {y, 1, 6}] ==\n",
        "  Table[Probability[X <= x, rule] Probability[Y <= y, rule],\n",
        "  {x, 1, 6}, {y, 1, 6}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "distU = DiscreteUniformDistribution[{1, 6}];\n",
        "distXY = TransformedDistribution[{Mod[U, 2], Mod[U, 3]},\n",
        "   Distributed[U, distU]];\n",
        "rule = Distributed[{X, Y}, distXY];\n",
        "Table[Probability[And[X <= x, Y <= y], rule], {x, 0, 1}, {y, 0, 2}] ==\n",
        " Table[Probability[X <= x, rule] Probability[Y <= y, rule],\n",
        "  {x, 0, 1}, {y, 0, 2}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "distX = BinomialDistribution[3, 1/2];\n",
        "distXY = TransformedDistribution[\n",
        "   {X, Piecewise[{{1, Or[X == 0, X == 3]}}, 2]}, Distributed[X, distX]];\n",
        "Covariance[distXY][[1, 2]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "rule = Distributed[{X, Y}, distXY];\n",
        "Table[Probability[And[X <= x, Y <= y], rule], {x, 0, 3}, {y, 1, 2}] ==\n",
        " Table[Probability[X <= x, rule] Probability[Y <= y, rule],\n",
        "  {x, 0, 3}, {y, 1, 2}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "c = Counts[Flatten[Table[{Max[x, y], Min[x, y]}, {x, 1, 6}, {y, 1, 6}], 1]]/36;\n",
        "dist = ProbabilityDistribution[Piecewise[KeyValueMap[{#2, {X, Y} == #1} &, c]],\n",
        "   {X, 1, 6, 1}, {Y, 1, 6, 1}];\n",
        "rule = Distributed[{X, Y}, dist];\n",
        "\n",
        "{Expectation[X + Y, rule],\n",
        " Expectation[X, rule] + Expectation[Y, rule]} (* 平均1 *)\n",
        "\n",
        "{distX, distY} = Table[MarginalDistribution[dist, i], {i, 2}];\n",
        "distXplusY = TransformedDistribution[X + Y, rule];\n",
        "\n",
        "{Mean[distXplusY], Mean[distX] + Mean[distY]} (* 平均2 *)\n",
        "{Variance[distXplusY],\n",
        " Variance[distX] + Variance[distY] + 2 Covariance[dist][[1, 2]]} (* 分散 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Manipulate[\n",
        " distY = BinomialDistribution[n, p];\n",
        " mu = Mean[distY]; sigma = StandardDeviation[distY];\n",
        " distZ = NormalDistribution[mu, sigma];\n",
        " Show[DiscretePlot[PDF[distY][x], {x, 0, n}], Plot[PDF[distZ][x], {x, 0, n}]],\n",
        " {{n, 15}, 1, 40, 1}, {{p, 4/10}, 0, 1}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "distX = UniformDistribution[]; distZ = NormalDistribution[];\n",
        "data = Table[Total[RandomVariate[distX, 12]] - 6, {10000}];\n",
        "Show[Histogram[data, Automatic, \"PDF\"], Plot[PDF[distZ][x], {x, -4, 4}]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist1 = NormalDistribution[0, 2]; dist2 = NormalDistribution[1, 1];\n",
        "TransformedDistribution[{X1 + X2 + 2, X1 + 3 X2 + 3},\n",
        " {Distributed[X1, dist1], Distributed[X2, dist2]}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "mu = {3, 6}; Sigma = {{5, 7}, {7, 13}};\n",
        "dist = MultinormalDistribution[mu, Sigma];\n",
        "Plot3D[PDF[dist][{Y1, Y2}], {Y1, -8, 14}, {Y2, -5, 17},\n",
        " PlotPoints -> 100, PlotRange -> All, AxesLabel -> Automatic]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "ContourPlot[PDF[dist][{Y1, Y2}], {Y1, -8, 14}, {Y2, -5, 17},\n",
        " PlotPoints -> 50, FrameLabel -> Automatic]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "sol = Solve[{Y1 == X1 + X2 + 2, Y2 == X1 + 3 X2 + 3}, {X1, X2}][[1]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "J = D[{X1, X2} /. sol, {{Y1, Y2}}];\n",
        "absj = Abs[Det[J]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f1 = PDF[dist1]; f2 = PDF[dist2];\n",
        "PDF[dist][{Y1, Y2}] == f1[X1] f2[X2] absj /. sol // Simplify"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "{MarginalDistribution[dist, 1], MarginalDistribution[dist, 2]}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = MultinormalDistribution[{u1, u2}, {{v1, 0}, {0, v2}}];\n",
        "d1 = MarginalDistribution[dist, 1]; d2 = MarginalDistribution[dist, 2];\n",
        "Simplify[CDF[dist][{x1, x2}] == CDF[d1][x1] CDF[d2][x2],\n",
        " And[v1 >= 0, v2 >= 0]]"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 10 推測統計",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = NormalDistribution[2, 3];\n",
        "data1 = Table[Mean[RandomVariate[dist, 5]], 10000];\n",
        "data2 = Table[Mean[RandomVariate[dist, 50]], 10000];\n",
        "{{Mean[data1], Variance[data1]}, {Mean[data2], Variance[data2]}}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Histogram[{data1, data2}, ChartLayout -> \"Row\"]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = NormalDistribution[2, 3];\n",
        "data1 = Table[Variance[RandomVariate[dist, 5]], 10000];\n",
        "data2 = Table[Variance[RandomVariate[dist, 50]], 10000];\n",
        "{{Mean[data1], Variance[data1]}, {Mean[data2], Variance[data2]}}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Histogram[{data1, data2}, ChartLayout -> \"Row\"]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "n = 4; mu = 5; sigma = 7; dist := NormalDistribution[mu, sigma];\n",
        "f[x_] := (n - 1) Variance[x]/sigma^2\n",
        "data = Table[f[RandomVariate[dist, n]], 10000];\n",
        "Show[Histogram[data, Automatic, \"PDF\"],\n",
        " Plot[PDF[ChiSquareDistribution[n - 1]][x], {x, 0, Max[data]}]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "n = 4; mu = 5; sigma = 7; ndist = NormalDistribution[mu, sigma];\n",
        "t = Function[{x}, (Mean[x] - mu)/Sqrt[Variance[x]/n]];\n",
        "data = Table[t[RandomVariate[ndist, n]], 10000];\n",
        "Show[Histogram[data, Automatic, \"PDF\"],\n",
        " Plot[PDF[StudentTDistribution[n - 1]][x], {x, -4.5, 4.5}]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "m = 5; muX = 2; sigmaX = 3; distX = NormalDistribution[muX, sigmaX];\n",
        "n = 7; muY = 3; sigmaY = 2; distY = NormalDistribution[muY, sigmaY];\n",
        "f[x_, y_] := (Variance[x]/sigmaX^2)/(Variance[y]/sigmaY^2)\n",
        "data = Table[f[RandomVariate[distX, m], RandomVariate[distY, n]], {10000}];\n",
        "Show[Histogram[data, Automatic, \"PDF\"],\n",
        " Plot[PDF[FRatioDistribution[m - 1, n - 1]][x], {x, 0, 7}]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[k, T];\n",
        "TransformedDistribution[T^2, Distributed[T, StudentTDistribution[k]]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "n = 15; p0 = 4/10; dist = BinomialDistribution[n, p0];\n",
        "tmp = Table[PDF[dist][x], {x, 0, n}];\n",
        "Total[Cases[tmp, p_ /; p <= PDF[dist][2]]] // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "CDF[dist][2] // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "n = 15; p0 = 4/10; dist = NormalDistribution[n p, Sqrt[n p (1 - p)]];\n",
        "2 CDF[dist /. p -> p0][2] // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "alpha = 5/100; InverseCDF[dist /. {p -> p0}, {alpha/2, 1 - alpha/2}] // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "N[Reduce[InverseCDF[dist, alpha/2] <= 2 <= InverseCDF[dist, 1 - alpha/2], p]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "pvalue[p0_] := With[{c = CDF[dist][2] /. p -> p0}, 2 Min[c, 1 - c]]\n",
        "Plot[pvalue[p0], {p0, 0, 1}]\n",
        "\n",
        "Plot[{InverseCDF[dist, alpha/2], InverseCDF[dist, 1 - alpha/2], 2},\n",
        " {p, 0, 1}, PlotStyle -> {Dashed, Thick, Dotted}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x = {24.2, 25.3, 26.2, 25.7, 24.4, 25.1, 25.6}; mu0 = 25;\n",
        "TTest[x, mu0]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "m = Mean[x]; s2 = Variance[x]; n = Length[x];\n",
        "t := (m - mu0)/Sqrt[s2/n];\n",
        "dist = StudentTDistribution[n - 1]; c = CDF[dist][t];\n",
        "2 Min[c, 1 - c]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "alpha = 5/100;\n",
        "{a, b} = InverseCDF[dist, {alpha/2, 1 - alpha/2}] // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Needs[\"HypothesisTesting`\"] (*「`」はシングルクォートではなくバッククォート *)\n",
        "MeanCI[x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[mu0]; Reduce[a <= t <= b, mu0]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "dist = StudentTDistribution[n - 1];\n",
        "Reduce[InverseCDF[dist, alpha/2] <= t <= InverseCDF[dist, 1 - alpha/2]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x = {25, 24, 25, 26}; y = {23, 18, 22, 28, 17, 25, 19, 16};\n",
        "TTest[{x, y}, 0, AlternativeHypothesis -> \"Greater\",\n",
        " VerifyTestAssumptions -> \"EqualVariance\" -> False]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "alpha = 5/100;\n",
        "m = Length[x]; n = Length[y]; sx2 = Variance[x]; sy2 = Variance[y];\n",
        "s2 = ((m - 1) sx2 + (n - 1) sy2)/(m + n - 2);\n",
        "T = (Mean[x] - Mean[y] - d)/Sqrt[s2 (1/m + 1/n)]; (* t統計量 *)\n",
        "t := T /. d -> 0                                  (* t値 *)\n",
        "df = m + n - 2;                                   (* 自由度 *)\n",
        "dist := StudentTDistribution[df];                 (* t分布 *)\n",
        "P := 1 - CDF[dist][t];                            (* P値 *)\n",
        "a := InverseCDF[dist, 1 - alpha];                 (* 採択域の上限 *)\n",
        "interval := Reduce[T <= a, d]                     (* 信頼区間 *)\n",
        "{t, P, a, interval} // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "T = (Mean[x] - Mean[y] - d)/Sqrt[sx2/m + sy2/n];\n",
        "df = (sx2/m + sy2/n)^2/((sx2/m)^2/(m - 1) + (sy2/n)^2/(n - 1)) // N;\n",
        "{t, P, a, interval} // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x = {25, 24, 25, 26}; y = {23, 18, 22, 28, 17, 25, 19, 16};\n",
        "VarianceTest[{x, y}, 1, \"HypothesisTestData\"][\"TestDataTable\"]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "m = Length[x]; n = Length[y]; dist = FRatioDistribution[m - 1, n - 1];\n",
        "F = Variance[x]/Variance[y]/r; f = F /. r -> 1;\n",
        "c = CDF[dist][f];\n",
        "{f, 2 Min[c, 1 - c]} // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "alpha = 5/100;\n",
        "{a, b} = InverseCDF[dist, {alpha/2, 1 - alpha/2}] // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Needs[\"HypothesisTesting`\"] (*「`」はシングルクォートではなくバッククォート *)\n",
        "VarianceRatioCI[x, y]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[a <= F <= b, r]"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 11 線形回帰分析",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "data = {{1, 2, 3}, {1, 3, 6}, {2, 5, 3}, {3, 7, 6}};\n",
        "model = LinearModelFit[data, {X1, X2}, {X1, X2}]\n",
        "model[\"BestFitParameters\"]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "model[1.5, 4]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x1 = {1, 3, 6, 10}; y = {7, 1, 6, 14};\n",
        "e = y - (b0 + b1 x1);\n",
        "L = e . e; (* 内積 *)\n",
        "FindMinimum[L, {{b0, 0}, {b1, 0}}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Minimize[L, {b0, b1}] (* 解析的な結果 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "L = Total[Abs[e]]; (* 差の絶対値の和 *)\n",
        "Minimize[L, {b0, b1}] // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "e = x1 - (y - b0)/b1;\n",
        "L = e . e;\n",
        "Minimize[L, {b0, b1}] // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "line = Module[{x1, y}, ImplicitRegion[y == b0 + b1 x1, {x1, y}]];\n",
        "L = Sum[RegionDistance[line, p]^2, {p, Thread[{x1, y}]}];\n",
        "Minimize[L, {b0, b1}] // Simplify // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "data = {{1, 2, 3}, {1, 3, 6}, {2, 5, 3}, {3, 7, 6}};\n",
        "X = DesignMatrix[data, {X1, X2}, {X1, X2}];\n",
        "y = data[[All, -1]];\n",
        "Inverse[Transpose[X] . X] . Transpose[X] . y"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "PseudoInverse[X] . y"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "b = {b0, b1, b2};\n",
        "L = (y - X . b) . (y - X . b);\n",
        "Reduce[{D[L, {b}] == 0 b}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "D[L, {b}] == -2 Transpose[X] . y + 2 Transpose[X] . X . b // Simplify"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "data = {{1, 2, 3}, {1, 3, 6}, {2, 5, 3}, {3, 7, 6}};\n",
        "model = LinearModelFit[data, {X1, X2}, {X1, X2}];\n",
        "model[\"RSquared\"]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "model[\"AdjustedRSquared\"]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x1 = {1, 3, 6, 10}; y = {7, 1, 6, 14}; data = Thread[{x1, y}];\n",
        "X = DesignMatrix[data, X1, X1];\n",
        "yh = X . PseudoInverse[X] . y;\n",
        "eh = y - yh; fh = yh - Mean[y]; g = y - Mean[y];\n",
        "R2 = 1 - eh . eh/g . g; N[R2]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "{Mean[eh] == 0,                                (* 特徴1 *)\n",
        " Mean[yh] == Mean[y],                          (* 特徴2 *)\n",
        " g . g == fh . fh + eh . eh,                   (* 特徴3 *)\n",
        " R2 == fh . fh/g . g,                          (* 特徴4 *)\n",
        " R2 == Correlation[y, yh]^2,                   (* 特徴5 *)\n",
        " 0 <= R2 <= 1,                                 (* 特徴6 *)\n",
        " Correlation[y, yh]^2 == Correlation[y, x1]^2} (* 特徴7 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "data = {{1, 2, 3}, {1, 3, 6}, {2, 5, 3}, {3, 7, 6}};\n",
        "n := Length[data]                                  (* サンプルサイズ *)\n",
        "p := Length[data[[1]]]                             (* 変数の個数 *)\n",
        "vars := Table[Subscript[x, i], {i, p - 1}]         (* 入力変数（記号）*)\n",
        "X := DesignMatrix[data, vars, vars]                (* 計画行列 *)\n",
        "y := data[[All, -1]]                               (* 出力変数の実現値 *)\n",
        "beta := Table[Subscript[\\[Beta], i - 1], {i, p}]   (* 回帰係数 *)\n",
        "epsilon := Table[Subscript[\\[Epsilon], i], {i, n}] (* 誤差項 *)\n",
        "Y := X . beta + epsilon                            (* 出力変数（確率変数）*)\n",
        "betah := PseudoInverse[X] . Y                      (* 回帰係数の推定量 *)\n",
        "betah // Simplify"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[sigma];\n",
        "udist = UniformDistribution[{-Sqrt[3] sigma, Sqrt[3] sigma}];\n",
        "udists = Table[Distributed[v, udist], {v, epsilon}];\n",
        "Expectation[betah, udists]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "ndist = NormalDistribution[0, sigma];\n",
        "ndists = Table[Distributed[v, ndist], {v, epsilon}];\n",
        "TransformedDistribution[betah, ndists] ==\n",
        " MultinormalDistribution[beta, sigma^2 Inverse[Transpose[X] . X]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "model := LinearModelFit[data, vars, vars];\n",
        "model[\"EstimatedVariance\"]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "e := Y - X . betah; RSS := e . e; s2 := RSS/(n - p)\n",
        "s2 // Simplify"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Expectation[s2, udists]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "tmp = Block[{sigma = 2},\n",
        "   dist = TransformedDistribution[Simplify[(n - p)  s2/sigma^2], ndists];\n",
        "   RandomVariate[dist, 10000]];\n",
        "cdist = ChiSquareDistribution[n - p];\n",
        "Show[Histogram[tmp, Automatic, \"PDF\"],\n",
        " Plot[PDF[cdist][x], {x, 0, 5}, PlotRange -> {0, 2}]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "uh := Transpose[A] . betah\n",
        "M := Transpose[A] . Inverse[Transpose[X] . X] . A\n",
        "r := MatrixRank[A]\n",
        "F := (uh - u) . Inverse[M] . (uh - u)/r/s2\n",
        "fdist := FRatioDistribution[r, n - p]\n",
        "pvalue := 1 - CDF[fdist, F]\n",
        "\n",
        "Y := y (* この先，実現値のみを扱う．*)\n",
        "A = Transpose[{{0, 1, 0}, {0, 0, 1}}]; u = {0, 0};\n",
        "{F, pvalue} // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "model[\"ParameterTable\"]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "u = {0};\n",
        "A = Transpose[{{1, 0, 0}}]; pvalue // N (* k = 0 *)\n",
        "A = Transpose[{{0, 1, 0}}]; pvalue // N (* k = 1 *)\n",
        "A = Transpose[{{0, 0, 1}}]; pvalue // N (* k = 2 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "s := Sqrt[s2 Diagonal[Inverse[Transpose[X] . X]]]\n",
        "s // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "t := betah/s\n",
        "t // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "tdist := StudentTDistribution[n - p]\n",
        "Table[2 Min[CDF[tdist][v], 1 - CDF[tdist][v]], {v, t}] // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "data = Transpose[{{35, 45, 55, 65, 75}, {114, 124, 143, 158, 166}}];\n",
        "alpha = 5/100; level := ConfidenceLevel -> 1 - alpha\n",
        "model[\"ParameterConfidenceIntervalTable\", level]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "tmp = InverseCDF[tdist, 1 - alpha/2];\n",
        "{betah - s tmp, betah + s tmp} // Transpose // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "cond := F <= InverseCDF[fdist, 1 - alpha]\n",
        "confint := Reduce[cond]\n",
        "A = Transpose[{{1, 0}}]; u = {beta0}; confint // N (* k = 0 *)\n",
        "A = Transpose[{{0, 1}}]; u = {beta1}; confint // N (* k = 1 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "tmp = model[\"ParameterConfidenceIntervals\", level];\n",
        "g1 = Graphics[{Gray, Apply[Rectangle, Transpose[tmp]]}];\n",
        "g2 = Graphics[model[\"ParameterConfidenceRegion\", level]];\n",
        "Show[g1, g2, AspectRatio -> 1, Frame -> True]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{1, 0}, {0, 1}}; u = {beta0, beta1};\n",
        "g3 = RegionPlot[ImplicitRegion[N[cond], Evaluate[u]]];\n",
        "Show[g1, g3, AspectRatio -> 1, Frame -> True]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "data = Transpose[{{35, 45, 55, 65, 75}, {114, 124, 143, 158, 166}}];\n",
        "g = Show[ListPlot[data], Plot[model[x1], {x1, 35, 75}],\n",
        "  Plot[Evaluate[model[\"MeanPredictionBands\", level]],\n",
        "   Evaluate[{vars[[1]], 35, 75}], PlotStyle -> Dashed]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Show[ListPlot[data], Plot[model[x1], {x1, 35, 75}],\n",
        "   Plot[Evaluate[model[\"SinglePredictionBands\", level]],\n",
        "     Evaluate[{vars[[1]], 35, 75}], PlotStyle -> Dashed]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{1}, {vars[[1]]}}; u = {Yp};\n",
        "Show[g, RegionPlot[Evaluate[cond],\n",
        "  Evaluate[{vars[[1]], 35, 75}], {Yp, 0, 200}, BoundaryStyle -> None]]"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 12 関数の極限と連続性",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[x_] := 2 x - 3\n",
        "Limit[f[x], x -> 1]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Limit[2 x - 3, x -> 1]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[x_] := Piecewise[{{x^2, x != 2}, {3, x == 2}}]\n",
        "Limit[f[x], x -> 2]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "g[x_] := (x^2 - 2)/(x - Sqrt[2])\n",
        "Limit[g[x], x -> Sqrt[2]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A := ForAll[epsilon, epsilon > 0, Exists[delta, delta > 0, B]];\n",
        "B := ForAll[x, Element[x, Reals],\n",
        "  Implies[0 < Norm[x - a] < delta, Norm[f[x] - alpha] < epsilon]]\n",
        "\n",
        "f[x_] := 2 x - 3; a = 1; alpha = -1;\n",
        "Reduce[A, Reals]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Simplify[Reduce[B, Reals], epsilon > 0]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[alpha];\n",
        "Reduce[A, Reals]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Limit[(1 + 1/x)^x, x -> Infinity]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Limit[1/x^2, x -> 0]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "{Limit[RealAbs[x]/x, x -> 0, Direction -> \"FromAbove\"],\n",
        " Limit[RealAbs[x]/x, x -> 0, Direction -> \"FromBelow\"]}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[f, g, x];\n",
        "f[x_] := Piecewise[{{RealAbs[x]/x, x != 0}}, Undefined]\n",
        "g[x_] := Piecewise[{{(x^2 - 2)/(x - Sqrt[2]), x != Sqrt[2]}}, Undefined]\n",
        "\n",
        "ResourceFunction[\"EnhancedPlot\"][f[x], {x, -1, 1}, \"FindExceptions\" -> True]\n",
        "ResourceFunction[\"EnhancedPlot\"][g[x], {x, 0, 2},  \"FindExceptions\" -> True]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "FunctionContinuous[{f[x], x != 0}, x]\n",
        "FunctionContinuous[{g[x], x != Sqrt[2]}, x]"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 13 微分",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[x_] := x^3\n",
        "f'[1]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "a = 1;\n",
        "Limit[(f[a + h] - f[a])/h, h -> 0]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[x_] := x^3\n",
        "f'[x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[x_] := x^3\n",
        "f1 = f'          (* 方法1 *)\n",
        "Derivative[1][f] (* 方法2 *)\n",
        "f2 = f'[x]       (* 方法1 *)\n",
        "D[f[x], x]       (* 方法2 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "{f1[1], f2 /. x -> 1}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "D[x^3, {x, 2}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[a, b, f, g];\n",
        "f[t_] := t^2\n",
        "g[x_] := a x + b\n",
        "Composition[f, g]'[x]                (* ① *)\n",
        "D[f[g[x]], x]                        (* ② *)\n",
        "(D[f[t], t] /. t -> g[x]) D[g[x], x] (* ③ *)\n",
        "f'[g[x]] g'[x]                       (* ④ *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "tmp = Series[Sin[x], {x, 0, 5}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Plot[Evaluate[{Normal[tmp], Sin[x]}], {x, -Pi, Pi}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "a = 0; Sum[Derivative[k][Sin][a] (x - a)^k/k!, {k, 0, 5}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[x_] := Sqrt[1 + x]\n",
        "Series[f[x], {x, 0, 7}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "fn[n_, a_, x_] := Sum[Derivative[k][f][a] (x - a)^k/k!, {k, 0, n},\n",
        "  GenerateConditions -> True]\n",
        "Reduce[f[x] == fn[Infinity, 0, x], Reals]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[x_] := Piecewise[{{Exp[-1/x^2], x != 0}}, 0]\n",
        "Plot[f[x], {x, -1, 1}]\n",
        "Reduce[f[x] == fn[Infinity, 0, x], Reals]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Derivative[k][f][0]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[x_] := x^3 - 12 x\n",
        "ResourceFunction[\"LocalExtrema\"][f[x], x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "sol = SolveValues[f'[x] == 0, x]\n",
        "Series[f[x], {x, sol[[1]], 2}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[a, delta];\n",
        "f[x_] := Piecewise[{{Exp[-1/x^2], x != 0}}, 0]\n",
        "Reduce[Exists[delta, delta > 0, ForAll[x, Element[x, Reals],\n",
        "   Implies[0 < Norm[x - a] < delta, f[a] < f[x]]]], Reals]"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 14 積分",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Integrate[-x^2 + 4 x + 1, {x, 1, 4}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[x_] := -x^2 + 4 x + 1\n",
        "Clear[x]; a = 1; b = 4; h = (b - a)/n;\n",
        "s = Sum[f[a + k h] h, {k, 1, n}] // Expand\n",
        "Limit[s, n -> Infinity]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Integrate[-t^2 + 4 t + 1, {t, a, x}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Integrate[-x^2 + 4 x + 1, x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[x, y];\n",
        "DSolveValue[y'[x] == -x^2 + 4 x + 1, y[x], x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "DSolveValue[{y'[x] == -x^2 + 4 x + 1, y[0] == 1}, y[x], x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "tmp = DSolveValue[y'[x] == -x y[x], y[x], x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[Integrate[tmp, {x, -Infinity, Infinity}] == 1]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[a, f, t, x];\n",
        "Function[x, Evaluate[Integrate[f[t], {t, a, x}]]]'"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "D[Integrate[f[t], {t, a, x}], x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "F = Integrate[-x^2 + 4 x + 1, x];\n",
        "(F /. x -> 4) - (F /. x -> 1)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Integrate[Log[Sin[x]], {x, 0, Pi/2}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[x_] := 1/(2 + Cos[x])\n",
        "F1x = Integrate[f[x], x];\n",
        "(F1x /. x -> 2 Pi) - (F1x /. x -> 0) (* 不正解 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "F2x = Integrate[f[t], {t, 0, x}, GenerateConditions -> True];\n",
        "(F2x /. x -> 2 Pi) - (F2x /. x -> 0) (* 正解 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "GraphicsRow[{Plot[F1x, {x, 0, 2 Pi}], Plot[F2x, {x, 0, 2 Pi}]}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Integrate[(p x + q)^100, x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "tmp = IntegrateChangeVariables[\n",
        "  Inactive[Integrate][(p x + q)^100, x], u, u == p x + q]\n",
        "Activate[tmp] /. u -> p x + q"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[x, y];\n",
        "IntegrateChangeVariables[\n",
        " Inactive[Integrate][1, {x, 0, t}], y, x == Sqrt[y]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Integrate[1/x^a, {x, 0, 1}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Integrate[1/x^a, {x, 1, Infinity}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Integrate[Exp[-x^2], {x, -Infinity, Infinity}]"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 15 多変数関数の微分積分",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x = {x1, x2}; f[{x1_, x2_}] := x1 x2^2/(x1^2 + x2^2)\n",
        "Limit[f[x], x -> {0, 0}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A := ForAll[epsilon, epsilon > 0, Exists[delta, delta > 0, B]];\n",
        "B := ForAll[Evaluate[x], Element[x, Reals],\n",
        "  Implies[0 < Norm[x - a] < delta, Norm[f[x] - alpha] < epsilon]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "a = {0, 0}; alpha = 0;\n",
        "Reduce[A, Reals]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[alpha];\n",
        "Reduce[A, Reals]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[x, y]; f[x_, y_] := x^2 y/(x^4 + y^2)\n",
        "Limit[f[x, y], {x, y} -> {0, 0}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[x, y, r, theta];\n",
        "{Limit[Limit[f[x, y], x -> 0], y -> 0],        (* ① *)\n",
        " Limit[Limit[f[x, y], y -> 0], x -> 0],        (* ② *)\n",
        " Limit[f[r Cos[theta], r Sin[theta]], r -> 0], (* ③ *)\n",
        " Limit[f[x, x^2], x -> 0]}                     (* ④ *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[{x1_, x2_}] := Piecewise[{{0, x1 == x2 == 0}}, x1 x2^2/(x1^2 + x2^2)]\n",
        "x = {x1, x2};\n",
        "FunctionContinuous[f[x], x]           (* 方法1 *)\n",
        "Limit[f[x], x -> {0, 0}] == f[{0, 0}] (* 方法2 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[x_, y_] := Piecewise[{{0, x == y == 0}}, x^2 y/(x^4 + y^2)]\n",
        "Clear[x, y]; FunctionContinuous[f[x, y], {x, y}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[x_, y_] := 2 - x^2 - y^2\n",
        "{D[f[x, y], x], D[f[x, y], y]}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[x_, y_] := 2 - x^2 - y^2\n",
        "{Derivative[1, 0][f], Derivative[0, 1][f]}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "g[{x1_, x2_}] := 2 - x1^2 - x2^2\n",
        "{Derivative[{1, 0}][g], Derivative[{0, 1}][g]}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "D[f[x, y], {{x, y}}]  (* 方法1 *)\n",
        "Grad[f[x, y], {x, y}] (* 方法2 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[x_, y_] := 2 x^3 + 5 x y + 2 y^2\n",
        "D[f[x, y], {{x, y}, 2}] // MatrixForm"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[f, F];\n",
        "f[{x1_, x2_}] := Sqrt[x1^2 + x2^2]\n",
        "x = {x1, x2}; a = {1, 1}; h = x - a;\n",
        "F[t_] := f[a + t h]\n",
        "expr := Normal[Series[F[t], {t, 0, 2}]] /. t -> 1\n",
        "expr // Simplify"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Block[{h = {h1, h2}}, expr /. Thread[h -> Map[HoldForm, x - a]]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "gradf = D[f[x], {x}] /. Thread[x -> a];\n",
        "H = D[f[x], {x, 2}] /. Thread[x -> a];\n",
        "f[a] + gradf . (x - a) + (x - a) . H . (x - a)/2 // Simplify"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x = {x1, x2}; f[{x1_, x2_}] := 2 x1^3 + x1 x2^2 + 5 x1^2 + x2^2\n",
        "ResourceFunction[\"LocalExtrema\"][f[x], x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "points := Solve[D[f[x], {x}] == 0 x, x, Reals]; (* 停留点 *)\n",
        "H := D[f[x], {x, 2}];                           (* ヘッセ行列 *)\n",
        "Table[With[{h = H /. p},                        (* 停留点でのヘッセ行列 *)\n",
        "  {p, f[x] /. p, Which[\n",
        "    PositiveDefiniteMatrixQ[h], -1,             (* 極小 *)\n",
        "    NegativeDefiniteMatrixQ[h], 1,              (* 極大 *)\n",
        "    IndefiniteMatrixQ[h], 0,                    (* 極値ではない *)\n",
        "    True, Null]}],                              (* 不明 *)\n",
        " {p, points}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x = {x1, x2}; f[{x1_, x2_}] := x1^2 + x2^4\n",
        "PositiveDefiniteMatrixQ[H /. Thread[x -> {0, 0}]] (* False *)\n",
        "ResourceFunction[\"LocalExtrema\"][f[x], x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "d = ImplicitRegion[And[0 <= x <= 1, 0 <= y <= x], {x, y}];\n",
        "f[x_, y_] := x^2 + y^2\n",
        "Integrate[f[x, y], Element[{x, y}, d]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Integrate[Integrate[f[x, y], {y, 0, x}], {x, 0, 1}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Integrate[Integrate[f[x, y], {x, y, 1}], {y, 0, 1}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[u, v, x, y];\n",
        "lhs = Inactive[Integrate][f[x, y], Element[{x, y}, d]]\n",
        "rhs = IntegrateChangeVariables[lhs, {u, v}, {x == 2 u, y == 3 v}]\n",
        "{Activate[lhs], Activate[rhs]}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "f[x_, y_] := x^2 + y^2\n",
        "{x, y} = {2 u, 3 v};\n",
        "J = D[{x, y}, {{u, v}}];\n",
        "detJ = Det[J]\n",
        "Integrate[Integrate[f[x, y] Abs[detJ], {v, 0, 2 u/3}], {u, 0, 1/2}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "{x, y} = {r Cos[theta], r Sin[theta]};\n",
        "J = D[{x, y}, {{r, theta}}];\n",
        "Det[J] // Simplify"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[x, y];\n",
        "lhs = Inactive[Integrate][Exp[-(x^2 + y^2)],\n",
        "  {y, -Infinity, Infinity}, {x, -Infinity, Infinity}]\n",
        "rhs = IntegrateChangeVariables[lhs, {r, theta}, \"Cartesian\" -> \"Polar\"]\n",
        "{Activate[lhs], Activate[rhs]}"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 16 ベクトル",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "a = {1/10 + 2/10, 1/10 + 2/10 - 3/10}; b = {3/10, 0};\n",
        "a == b"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "100 {1, 2} + 10 {3, 1}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "a = {3, 4};\n",
        "Norm[a]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[x, y]; a = {x, y}; Sqrt[a . a]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Simplify[Norm[{x, y}], Element[x | y, Reals]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "a = {3, 4};\n",
        "Normalize[a]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "a = {1, 0}; b = {1, 1};\n",
        "ArcCos[a . b/(Norm[a] Norm[b])]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "VectorAngle[a, b]"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 17 行列",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "MatrixForm[A = {{1, 2, 0}, {0, 3, 4}}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x = {5, 7}; DiagonalMatrix[x] // MatrixForm"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "SymmetricMatrixQ[{{1, 2}, {2, 3}}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "MatrixForm[A = {{11, 12, 13}, {21, 22, 23}, {31, 32, 33}}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A[[1 ;; 2, 1 ;; 2]] // MatrixForm"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A[[All, 3]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A[[All, {3}]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A[[2, All]] (* 方法1 *)\n",
        "A[[2]]      (* 方法2 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A[[{2}, All]] (* 方法1 *)\n",
        "A[[{2}]]      (* 方法2 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "10 {{2, 3}, {5, 7}}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "{{10, 20}, {30, 40}} + {{2, 3}, {4, 5}}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{2, 3}, {5, 7}}; B = {{1, 2}, {3, 4}};\n",
        "A . B"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{2, 3}, {5, 7}}; B = {{1, 2, 3}, {4, 5, 6}}; S = A . B;\n",
        "{p, q} = Dimensions[A]; {r, s} = Dimensions[B];\n",
        "S1 = Table[Table[A[[i, All]] . B[[All, j]], {j, 1, s}], {i, 1, p}]; (* ① *)\n",
        "S2 = Sum[A[[All, {j}]] . B[[{j}, All]], {j, 1, q}];                 (* ② *)\n",
        "S3 = Transpose[Table[A . b, {b, Transpose[B]}]];                    (* ③ *)\n",
        "S4 = Table[a . B, {a, A}];                                          (* ④ *)\n",
        "{S == S1, S == S2, S == S3, S == S4}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[a1, a2, x1, x2, p, q, r, s];\n",
        "x = {x1, x2}; a = {a1, a2};\n",
        "G = {{p, q}, {q, s}}; A = {{p, q}, {r, s}};\n",
        "D[a . x, {x}] == a\n",
        "D[x . G . x, {x}] == 2 G . x // Simplify\n",
        "D[(A . x) . (A . x), {x}] == 2 Transpose[A] . A . x // Simplify"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Det[{{3, 2}, {1, 2}}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "RegionMeasure[Parallelepiped[{0, 0}, {{3, 1}, {2, 2}}]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "RegionMeasure[Parallelepiped[{0, 0, 0}, {{2, 1, 0}, {0, 2, 1}, {1, 1, 1}}]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Inverse[{{2, 3}, {5, 7}}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{3, 2}, {1, 2}}; b = {8, 4};\n",
        "Inverse[A] . b"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "RowReduce[{{4, 2, 8}, {2, 1, 4}}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{2, 0, 2}, {0, 2, -2}, {2, 2, 0}};\n",
        "MatrixRank[A]"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 18 ベクトル空間",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "a1 = {3, 1}; a2 = {2, 2};\n",
        "ResourceFunction[\"LinearlyIndependent\"][{a1, a2}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[c1 a1 + c2 a2 == {0, 0}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{1, 0, 1}, {1, 1, 0}, {0, 1, -1}};\n",
        "ResourceFunction[\"ColumnSpace\"][A][\"Basis\"]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{1, 0, 1}, {1, 1, 0}, {0, 1, -1}};\n",
        "tmp = ResourceFunction[\"ColumnSpace\"][A];\n",
        "Qt = Orthogonalize[tmp[\"Basis\"]]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Q = Transpose[Qt];\n",
        "Qt . Q"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{1, 2}, {1, 2}, {0, 0}}; B = {{1, 0}, {1, 1}, {0, 1}};\n",
        "{tQa, Ra} = QRDecomposition[A]; Qa = Transpose[tQa]; (* 転置が必要 *)\n",
        "{tQb, Rb} = QRDecomposition[B]; Qb = Transpose[tQb]; (* 転置が必要 *)\n",
        "{MatrixForm[Qa], MatrixForm[Ra], A == Qa . Ra,\n",
        " MatrixForm[Qb], MatrixForm[Rb], B == Qb . Rb}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "qrd[A_] := Module[{m, n, u = Transpose[A], idx = {}, s, Q},\n",
        "  {m, n} = Dimensions[A];\n",
        "  Do[Do[u[[i]] = Simplify[u[[i]] - A[[All, i]] . u[[j]] u[[j]]], {j, 1, i - 1}];\n",
        "   s = Chop[Norm[u[[i]]]];\n",
        "   If[s != 0, u[[i]] /= s; AppendTo[idx, i]], {i, 1, n}];\n",
        "  Q = If[Length[idx] != 0, Transpose[u[[idx]]], IdentityMatrix[m]];\n",
        "  {Q, Transpose[Q] . A}]\n",
        "\n",
        "A = {{1, 2}, {1, 2}, {0, 0}}; B = {{1, 0}, {1, 1}, {0, 1}};\n",
        "Map[MatrixForm, qrd[A]] // Simplify (* 動作確認 *)\n",
        "Map[MatrixForm, qrd[B]] // Simplify (* 動作確認 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "B = {{1, 0}, {1, 1}, {0, 1}};\n",
        "{Q, R} = qrd[B];                              (* QR分解 *)\n",
        "tol = 10^-10;\n",
        "e = IdentityMatrix[Dimensions[Q][[2]]];\n",
        "{Chop[N[Transpose[Q] . Q] - e, tol] == 0 e,   (* ① *)\n",
        " UpperTriangularMatrixQ[R, Tolerance -> tol], (* ② *)\n",
        " Chop[N[B] - Q . R, tol] == 0 B}              (* ③ *)\n",
        "(* 誤った転置を検出できないから，①でOrthogonalMatrixQは使えない．*)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{1, 0}, {1, 1}, {0, 1}};\n",
        "NullSpace[Transpose[A]]\n",
        "NullSpace[Transpose[N[A]]] (* 正規直交基底 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{1, 0}, {1, 1}, {0, 1}};\n",
        "basis1 = Orthogonalize[Transpose[A]];            (* 列空間 *)\n",
        "basis2 = Orthogonalize[NullSpace[Transpose[A]]]; (* 直交補空間 *)\n",
        "MatrixForm[Q = Transpose[Join[basis1, basis2]]]\n",
        "Transpose[Q] . Q == IdentityMatrix[3]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{a, b}, {c, d}};\n",
        "f[x_] := A . x\n",
        "R = ParametricRegion[{x, y}, {{x, s, s + u}, {y, t, t + u}}];\n",
        "Rp = TransformedRegion[R, f];\n",
        "{RegionMeasure[Rp], Abs[Det[A]] u^2}"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 19 固有値と固有ベクトル",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{5, 6, 3}, {0, 9, 2}, {0, 6, 8}}; (* 固有ベクトル（絶対値の降順） *)\n",
        "{vals, vecs} = Eigensystem[N[A]]       (* 近似値：固有ベクトル（正規） *)\n",
        "{vals, vecs} = Eigensystem[A]          (* 厳密値：固有ベクトル（非正規）*)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "V = Transpose[vecs]; A . V == V . DiagonalMatrix[vals]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{5, 6, 3}, {0, 9, 2}, {0, 6, 8}}; n = Length[A];\n",
        "SolveValues[Det[x IdentityMatrix[n] - A] == 0, x]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "NullSpace[5 IdentityMatrix[n] - A]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "S = {{2, 2, -2}, {2, 5, -4}, {-2, -4, 5}};\n",
        "{Q, L, V} = SingularValueDecomposition[S];\n",
        "{MatrixForm[Q], MatrixForm[L],\n",
        " S == Q . L . Transpose[Q] == V . L . Transpose[V]}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "S = {{2, 2, -2}, {2, 5, -4}, {-2, -4, 5}};\n",
        "{vals, vecs} = Eigensystem[S];           (* ① *)\n",
        "Q = Transpose[Orthogonalize[vecs]];      (* ②, ③ *)\n",
        "L = DiagonalMatrix[vals];                (* ④ *)\n",
        "Chop[N[S] - Q . L . Transpose[Q]] == 0 S (* 近似的な比較 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "PositiveSemidefiniteMatrixQ[{{4, 2}, {2, 1}}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{4, 2}, {2, 1}};\n",
        "AllTrue[Eigenvalues[A], NonNegative]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "x1 = {1, 3, 6, 10}; y = {7, 1, 6, 14}; X = Transpose[{x1, y}];\n",
        "n = Length[X]; M = ConstantArray[1, {n, n}]/n;\n",
        "A = X - M . X;\n",
        "MatrixForm[S = Transpose[A] . A]\n",
        "v = Eigenvectors[N[S], 1][[1]] (* 最大固有値に対応する固有ベクトル *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Reduce[Det[{v, {xp - Mean[x1], yp - Mean[y]}}] == 0, yp] // N"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "{U, L, V} = SingularValueDecomposition[A]; (* 特異値分解 *)\n",
        "V[[All, 1]] // N                           (* Vの第1列（求めるもの）*)\n",
        "s2 = Diagonal[L]^2;                        (* 特異値の2乗 *)\n",
        "Accumulate[s2]/Total[s2] // N              (* 累積寄与率（後述）*)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "X = N[Transpose[{{1, 3, 6, 10}, {7, 1, 6, 14}}]];\n",
        "t = Transpose;\n",
        "MatrixForm[P = PrincipalComponents[X]]                 (* 主成分スコア *)\n",
        "r = MatrixRank[P]; Pr = P[[All, ;; r]]; tPr = t[Pr];\n",
        "MatrixForm[tVr1 = Inverse[tPr . Pr] . tPr . X]         (* 主成分（方法1）*)\n",
        "MatrixForm[tVr2 = (PseudoInverse[P] . X)[[;; r, All]]] (* 主成分（方法2）*)\n",
        "tVr1[[1]]                        (* 第1主成分（求めるもの）（方法1）*)\n",
        "tVr2[[1]]                        (* 第1主成分（求めるもの）（方法2）*)\n",
        "s2 = Diagonal[Transpose[P] . P]; (* 特異値の2乗 *)\n",
        "Accumulate[s2]/Total[s2]         (* 累積寄与率（後述）*)"
      ]
    },
    {
      "cell_type": "markdown",
      "source": "# 20 特異値分解と擬似逆行列",
      "metadata": {}
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{1, 0}, {1, 1}, {0, 1}};\n",
        "{U, S, V} = SingularValueDecomposition[A]; tV = Transpose[V];\n",
        "{Map[MatrixForm, {U, S, tV}], A == U . S . tV}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "url = \"https://github.com/taroyabuki/comath/raw/main/images/boy.jpg\";\n",
        "A = ImageData[ColorConvert[Import[url], \"Grayscale\"]]; (* 画像の行列への変換 *)\n",
        "{U, S, V} = SingularValueDecomposition[A];             (* 特異値分解 *)\n",
        "k = 52;\n",
        "Ak = U[[All, ;; k]] . S[[;; k, ;; k]] . Transpose[V[[All, ;; k]]]; (* 近似 *)\n",
        "B = (Ak - Min[Ak])/(Max[Ak] - Min[Ak]); (* 数値を0～1にする．*)\n",
        "GraphicsRow[{Image[A], Image[B]}]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "nonzero[x_, tol_ : 10^-10] := Chop[x, tol] != 0\n",
        "svd2[A_] := Module[{diag = DiagonalMatrix, eye = IdentityMatrix, t = Transpose,\n",
        "   gs = Orthogonalize, m, n, G, vals, vecs, s, r, Sr, S, Vr, V, Ur, U},\n",
        "  {m, n} = Dimensions[A]; G = t[A] . A;                      (* ① *)\n",
        "  {vals, vecs} = Eigensystem[G];                             (* ② *)\n",
        "  s = Sqrt[Select[vals, nonzero]]; r = Length[s];            (* ③ *)\n",
        "  If[r != 0,\n",
        "   Sr = diag[s, 0, {r, r}];                                  (* ④ *)\n",
        "   Vr = t[gs[Take[vecs, r]]];                                (* ⑤ *)\n",
        "   Ur = A . Vr . diag[1/s, 0, {r, r}];                       (* ⑥ *)\n",
        "   S = diag[s, 0, {m, n}];                                   (* ⑦ *)\n",
        "   V = If[n == r, Vr, Join[Vr, t[gs[NullSpace[t[Vr]]]], 2]]; (* ⑧ *)\n",
        "   U = If[m == r, Ur, Join[Ur, t[gs[NullSpace[t[Ur]]]], 2]], (* ⑨ *)\n",
        "   (* else *)\n",
        "   S = 0 A; V = eye[n]; U = eye[m];\n",
        "   Sr = {{0}}; Vr = V[[All, {1}]]; Ur = U[[All, {1}]]];\n",
        "  {Ur, Sr, Vr, U, S, V}]\n",
        "\n",
        "A = {{1, 0}, {1, 1}, {0, 1}}; Map[MatrixForm, svd2[A]] (* 動作確認 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "tol = 10^-10;\n",
        "isOrtho[A_] := With[{e = IdentityMatrix[Dimensions[A][[2]]]},\n",
        "  Chop[Transpose[A] . A - e, tol] == 0 e]\n",
        "isDiagDesc[A_] := With[{d = Diagonal[A]}, d == Sort[Abs[d], Greater]]\n",
        "t = Transpose;\n",
        "\n",
        "A = {{1, 0}, {1, 1}, {0, 1}};\n",
        "{Ur, Sr, Vr, U, S, V} = svd2[A];                   (* 特異値分解 *)\n",
        "{isOrtho[Ur], isOrtho[Vr], isOrtho[U], isOrtho[V], (* ① *)\n",
        " SquareMatrixQ[U], SquareMatrixQ[V],               (* ② *)\n",
        " isDiagDesc[Sr], isDiagDesc[S],                    (* ③ *)\n",
        " Chop[N[A] - Ur . Sr . t[Vr], tol] == 0 A,         (* ④-1 *)\n",
        " Chop[N[A] - U . S . t[V], tol] == 0 A}            (* ④-2 *)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "Clear[\"Global`*\"];"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{1, 0}, {1, 1}, {0, 1}}; PseudoInverse[A]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "A = {{1, 0}, {1, 1}, {0, 1}}; b = {2, 0, 2};\n",
        "PseudoInverse[A] . b"
      ]
    }
  ],
  "nbformat": 4,
  "nbformat_minor": 5,
  "metadata": {
    "kernelspec": {
      "name": "wolframlanguage14.0",
      "display_name": "Wolfram Language 14.0"
    }
  }
}