PCA analysis
2022-09-14
## -- Attaching packages --------------------------------------- tidyverse 1.3.1 --## v ggplot2 3.3.6 v purrr 0.3.4
## v tibble 3.1.7 v dplyr 1.0.9
## v tidyr 1.2.0 v stringr 1.4.0
## v readr 2.1.2 v forcats 0.5.1## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()##
## 次のパッケージを付け加えます: 'plotly'## 以下のオブジェクトは 'package:ggplot2' からマスクされています:
##
## last_plot## 以下のオブジェクトは 'package:stats' からマスクされています:
##
## filter## 以下のオブジェクトは 'package:graphics' からマスクされています:
##
## layout## 要求されたパッケージ xts をロード中です## 要求されたパッケージ zoo をロード中です##
## 次のパッケージを付け加えます: 'zoo'## 以下のオブジェクトは 'package:base' からマスクされています:
##
## as.Date, as.Date.numeric##
## 次のパッケージを付け加えます: 'xts'## 以下のオブジェクトは 'package:dplyr' からマスクされています:
##
## first, last## 要求されたパッケージ TTR をロード中です## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo##
## 次のパッケージを付け加えます: 'magrittr'## 以下のオブジェクトは 'package:purrr' からマスクされています:
##
## set_names## 以下のオブジェクトは 'package:tidyr' からマスクされています:
##
## extract##
## 次のパッケージを付け加えます: 'scales'## 以下のオブジェクトは 'package:purrr' からマスクされています:
##
## discard## 以下のオブジェクトは 'package:readr' からマスクされています:
##
## col_factor## Welcome! Want to learn more? See two factoextra-related books at https://goo.gl/ve3WBaJGB Analysis[dataset source:財務省]
jgb2 <- read_excel("jgbcm_all.xlsx")
head(jgb2)## # A tibble: 6 x 16
## `1` `1Y` `2Y` `3Y` `4Y` `5Y` `6Y` `7Y` `8Y` `9Y`
## <dttm> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1974-09-24 00:00:00 10.3 9.36 8.83 8.52 8.35 8.29 8.24 8.12 8.13
## 2 1974-09-25 00:00:00 10.3 9.36 8.83 8.52 8.35 8.29 8.24 8.12 8.13
## 3 1974-09-26 00:00:00 10.3 9.37 8.83 8.52 8.35 8.29 8.24 8.12 8.13
## 4 1974-09-27 00:00:00 10.3 9.37 8.83 8.52 8.35 8.29 8.24 8.12 8.13
## 5 1974-09-28 00:00:00 10.4 9.37 8.83 8.52 8.35 8.29 8.24 8.12 8.13
## 6 1974-09-30 00:00:00 10.4 9.37 8.84 8.52 8.35 8.29 8.24 8.12 8.13
## # ... with 6 more variables: `10Y` <dbl>, `15Y` <dbl>, `20Y` <dbl>,
## # `25Y` <dbl>, `30Y` <dbl>, `40Y` <dbl>jgb2[is.na(jgb2)] <- 0
jgb3 <- jgb2[,-1]
### PCA analysis[multivariate analysis]
### PCA, (M)CA,FAMD,MFA,HCPC,factoextraResults
jgb_pca <- princomp((jgb3), cor = TRUE)
summary(jgb_pca)## Importance of components:
## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5
## Standard deviation 3.1433166 1.5537198 1.06638985 0.77222681 0.66221399
## Proportion of Variance 0.6586959 0.1609363 0.07581249 0.03975562 0.02923516
## Cumulative Proportion 0.6586959 0.8196323 0.89544477 0.93520039 0.96443554
## Comp.6 Comp.7 Comp.8 Comp.9
## Standard deviation 0.49534812 0.40145367 0.256220464 0.194928399
## Proportion of Variance 0.01635798 0.01074434 0.004376595 0.002533139
## Cumulative Proportion 0.98079353 0.99153786 0.995914458 0.998447597
## Comp.10 Comp.11 Comp.12 Comp.13
## Standard deviation 0.126705339 0.0629395558 0.0457497667 2.727444e-02
## Proportion of Variance 0.001070283 0.0002640925 0.0001395361 4.959302e-05
## Cumulative Proportion 0.999517880 0.9997819726 0.9999215086 9.999711e-01
## Comp.14 Comp.15
## Standard deviation 1.559118e-02 1.379820e-02
## Proportion of Variance 1.620566e-05 1.269269e-05
## Cumulative Proportion 9.999873e-01 1.000000e+00#Note that the `echo = FALSE` parameter was added to the code chunk to prevent printing of the R code that generated the plot.biplot(jgb_pca, main = "Biplot", scale = 0)
### compute standard deviation, variance, variance explained for each
principal component
jgb_pca$sdev## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6 Comp.7
## 3.14331656 1.55371981 1.06638985 0.77222681 0.66221399 0.49534812 0.40145367
## Comp.8 Comp.9 Comp.10 Comp.11 Comp.12 Comp.13 Comp.14
## 0.25622046 0.19492840 0.12670534 0.06293956 0.04574977 0.02727444 0.01559118
## Comp.15
## 0.01379820jgb_pca.var <- jgb_pca$sdev^2
jgb_pca.var## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6
## 9.8804389838 2.4140452405 1.1371873041 0.5963342508 0.4385273645 0.2453697621
## Comp.7 Comp.8 Comp.9 Comp.10 Comp.11 Comp.12
## 0.1611650453 0.0656489261 0.0379970807 0.0160542429 0.0039613877 0.0020930411
## Comp.13 Comp.14 Comp.15
## 0.0007438953 0.0002430849 0.0001903903#proportion of variance for screen plot
propve <- jgb_pca.var/sum(jgb_pca.var)
#plot variance explained for each p.c.
plot(propve, xlab = "principal component",
ylab = "proportion of variance explained",
ylim = c(0, 1), type = "b",
main = "screen plot")
factor.loadings <- jgb_pca$loadings[,1:3]
legend.loadings <- c("パラレル", "ツイスト", "バタフライ")
par(xaxt="n")
matplot(factor.loadings,
xlab = "Term", ylab="Factor Loadings")
#error
#legend(4,max(factor.loadings),legend=legend.loadings,col=1:3,lty=1,lwd=3)
#par(xaxt="s")
#axis(1,1:length(label.term),label.term)