noteで公表。

マトリックスの内部の計算

a <- matrix(1:12,nrow=3, byrow = F)
a
##      [,1] [,2] [,3] [,4]
## [1,]    1    4    7   10
## [2,]    2    5    8   11
## [3,]    3    6    9   12
b <- a[2,] + a[3,]
b
## [1]  5 11 17 23
c <- rbind(a,b)
c
##   [,1] [,2] [,3] [,4]
##      1    4    7   10
##      2    5    8   11
##      3    6    9   12
## b    5   11   17   23
d <- c[c(1,4),]
d
##   [,1] [,2] [,3] [,4]
##      1    4    7   10
## b    5   11   17   23
d <- c[-c(2,3),]
d
##   [,1] [,2] [,3] [,4]
##      1    4    7   10
## b    5   11   17   23
d <- c[c(-2,-3),]
d
##   [,1] [,2] [,3] [,4]
##      1    4    7   10
## b    5   11   17   23
rownames(d)[1] <- "a(r1)"
d
##       [,1] [,2] [,3] [,4]
## a(r1)    1    4    7   10
## b        5   11   17   23

2×2表

mydata <- matrix(c(70,30,45,5),2,byrow = T)
mydata
##      [,1] [,2]
## [1,]   70   30
## [2,]   45    5
addmargins(mydata, margin = 1:2)
##            Sum
##      70 30 100
##      45  5  50
## Sum 115 35 150
par(family="HiraKakuProN-W3")
barplot(mydata,xlab="性別",names=c("男","女"),col=c("blue","red"))

par(family="HiraKakuProN-W3")
barplot(mydata,names=c("男","女"),col=c("blue","red"), beside=T)

# 行パーセントの計算
m1 <- prop.table(mydata,1)
rownames(m1) <- c("male","female")
colnames(m1) <- c("yes","no")
m1
##        yes  no
## male   0.7 0.3
## female 0.9 0.1
m2 <- t(m1)*100
m2
##     male female
## yes   70     90
## no    30     10
barplot(m2,xlab="sex",col=c("blue","red"),legend=T,args.legend = list(x = "bottomright"))

barplot(m2,horiz=T,xlab="respose",col=c("blue","red"),legend=T,args.legend = list(x = "bottomright"))

m2
##     male female
## yes   70     90
## no    30     10
m3 <- m2[,c(2,1)]
m3
##     female male
## yes     90   70
## no      10   30
barplot(m3,horiz=T,xlab="respose",col=c("blue","red"),legend=T,args.legend = list(x = "bottomright"))

m_ <-  matrix(c(100-0.12,100-0.27,0.12,0.27),2,byrow = T)
barplot(m_,horiz=T,xlab="respose",col=c("blue","red"),legend=T,args.legend = list(x = "bottomright"))

p1 <- 0.12
p2 <- 0.27
p <-c(p1,p2)
barplot(p,ylim= c(0,0.3))

m_ <-  matrix(c(336,60306-336,24,6816-24),2,byrow = T)
m_ <- 100*prop.table(m_,1)
m_
##           [,1]     [,2]
## [1,] 0.5571585 99.44284
## [2,] 0.3521127 99.64789
rownames(m_) <- c("非着用","着用")
colnames(m_) <- c("致死","負傷")
par(family= "HiraKakuProN-W3")
barplot(t(m_),horiz = T, xlab = "致死率",ylab = c("リスク要因の有無"))

p1 <- 0.35
p2 <- 0.56
p <-c(p1,p2)
par(family= "HiraKakuProN-W3")
barplot(p,ylim= c(0,0.6) ,ylab = "致死率(単位：パーセント）",xlab=c("左：着用","右：非着用"))

mosaicplot(m2,color = T,main = "")

mosaicplot(t(m2),color = T,main = "")

m2
##     male female
## yes   70     90
## no    30     10
chisq.test(m2)
## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  m2
## X-squared = 11.281, df = 1, p-value = 0.0007829
chisq.test(m2, simulate.p.value=TRUE)
## 
##  Pearson's Chi-squared test with simulated p-value (based on 2000
##  replicates)
## 
## data:  m2
## X-squared = 12.5, df = NA, p-value = 0.0009995
fisher.test(m2)
## 
##  Fisher's Exact Test for Count Data
## 
## data:  m2
## p-value = 0.0006504
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
##  0.1063153 0.5936350
## sample estimates:
## odds ratio 
##  0.2609721

Everitt

x <- matrix(seq(12),nrow=4)
dimnames(x) <- list(c("R1","R2","R3","R4"), c("C1","C2","C3"))
x
##    C1 C2 C3
## R1  1  5  9
## R2  2  6 10
## R3  3  7 11
## R4  4  8 12
x <- matrix(seq(12),nrow=4)
dimnames(x) <- list(paste("row", seq(4)), paste("col", seq(3)))
x
##       col 1 col 2 col 3
## row 1     1     5     9
## row 2     2     6    10
## row 3     3     7    11
## row 4     4     8    12
height <- c(50,70,45,80,100)
weight <- c(120,140,100,200,190)
age <- c(20,40,41,31,33)
names <- c("Bob","Ted","Alice","Mary","Sue")
sex <- c("Male","Male","Female","Female","Female")
data <- data.frame(names,sex,height,weight,age)
data
data$age
## [1] 20 40 41 31 33
data[,"age"]
## [1] 20 40 41 31 33
data[,5]
## [1] 20 40 41 31 33
data[3:5,5]
## [1] 41 31 33
data[-c(1:2),5]
## [1] 41 31 33
data[c(-1,-2),5]
## [1] 41 31 33
table(data$sex)
## 
## Female   Male 
##      3      2
table(data$height, data$sex)
##      
##       Female Male
##   45       1    0
##   50       0    1
##   70       0    1
##   80       1    0
##   100      1    0
breaks <- seq(40,100,by=20)
result <- table(cut(data$height,breaks))
result
## 
##  (40,60]  (60,80] (80,100] 
##        2        2        1
pie(result, labels = c("A","B","C"), col = rainbow(3))

hist(data$height,breaks)

freq <- cut(data$height, breaks=breaks, 
                right=F, include.lowest=TRUE)
table(freq)
## freq
##  [40,60)  [60,80) [80,100] 
##        2        1        2
hist(data$height, breaks = breaks,right = F ,include.lowest = T)

散布図

plot(x=data$height, y=data$weight)
result <- lm(weight ~ height, data = data)
abline(result)

data
apply(data[,c(3,4,5)],2,mean) # 平均
## height weight    age 
##     69    150     33

計算例

data <- c(45,56,51,68,53,96,56,92,47,44,71,77,35,79,65,60,45,42,52,96,26,42,47,87,73,74,23,39,70,74,40,64)
data
##  [1] 45 56 51 68 53 96 56 92 47 44 71 77 35 79 65 60 45 42 52 96 26 42 47 87 73
## [26] 74 23 39 70 74 40 64
length(data)
## [1] 32
summary(data)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   23.00   44.75   56.00   59.03   73.25   96.00
sd(data)
## [1] 19.47784
breaks <- seq(0,100,by=20)
hist(data, breaks = breaks,right = F ,include.lowest = T)

breaks <- seq(0,100,by=5)
hist(data, breaks = breaks,right = F ,include.lowest = T)

boxplot(data)

ベクトルに用いる関数

x <- 1:100
sum(x)
## [1] 5050
summary(x)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    1.00   25.75   50.50   50.50   75.25  100.00
sd(x)
## [1] 29.01149
round(29.01149,digits=2)
## [1] 29.01
boxplot(x)

x <- rnorm(1000, mean=0, sd=1)
sum(x)
## [1] 24.50679
summary(x)
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
## -2.79099 -0.67781  0.01161  0.02451  0.72914  3.00760
round(sd(x),digits=2)
## [1] 0.99
hist(x)

boxplot(x)