abikoushi

library(readr)
library(dplyr)
library(ggplot2)

## data is available from
## https://www.data.jma.go.jp/risk/obsdl/index.php

dat = read_csv("./Downloads/data.csv", skip = 6, 
             col_names = c("date","Temp","hinshitsu","kinitsu"))

abikoushi

母比率の差の検定と分割表の独立性の検定が同じになること

次のように分割表が与えられたとき，

成功	失敗	合計
a	c	$n_1$
b	d	$n_2$
$m_1$	$m_2$	n

abikoushi

#佐藤俊哉『宇宙怪人しまりす　統計よりも重要なことを学ぶ』（朝倉書店）より
X = matrix(c(58, 22,
            62, 38), byrow = TRUE, nrow = 2, ncol = 2)
print(X)

res_chisq = chisq.test(X, correct = FALSE)
print(res_chisq$p.value)
#[1] 0.1375639

x = X[,1]

abikoushi

#include <chrono>
#include <iostream>
#include <fstream>
#include <mutex>

class TimerLogger {
public:
  TimerLogger(const std::string& label, const std::string& filename = "timelog.csv")
    : label_(label), filename_(filename) {
    start_ = std::chrono::high_resolution_clock::now();

abikoushi

data {
  int<lower=0> N;
  int<lower=0> D;
  int<lower=0> R;
  matrix[N,D] Y;
}
parameters {
  row_vector[D] mu;
  matrix[D,R] W;
  real<lower=0> sig2;

abikoushi

using Distributions
using CairoMakie
#https://docs.makie.org/stable/

xv = 0:0.1:20
f = Figure()
ax = Axis(f[1, 1])
l1 = lines!(ax, xv, pdf.(Gumbel(0,1), xv))
l2 = lines!(ax, xv, pdf.(Gumbel(1,1), xv), linestyle=:dash)
l3 = lines!(ax, xv, pdf.(Gumbel(0,2), xv), linestyle=:dot)

abikoushi

using Random
using Distributions
using CairoMakie
using FileIO

beta = [1, 2]

rng = Random.default_rng(1234)
X = [ones(10) randn(rng, 10)]
y = X * beta + randn(rng, 10)

abikoushi

makeLmat <- function (h){
    x = seq(0,1,by=h)
    N=length(x)
    invh2 <- 1/(h^2)
    res = diag(invh2, N)
    res[cbind(1:(N-1),2:N)] <- -0.5*invh2
    res[cbind(2:N,1:(N-1))] <- -0.5*invh2
    return(list(L=res,x=x))
}

abikoushi

library(ggplot2)
library(dplyr)
library(mvtnorm)
set.seed(1)
x1 = rnorm(10, 0, 1)
x2 = x1+rnorm(10, 0, 0.1)
y = x1+x2 + rnorm(10, 0, 1)

print(cor(x1,x2))

abikoushi

library(deSolve)

SEIRmod <- function(Time, State, Pars) {
  with(as.list(c(State, Pars)), {
    dS <- - beta*I*S
    dE <- beta*I*S - alpha*E
    dI <- alpha*E - gamma*I
    dR <- gamma*I
    return(list(c(dS, dE, dI, dR)))
  })