library(MASS)
capas <- array(NA, dim=c(50,2,99))
for (r in 1:99) {
z <- mvrnorm(50,mu=c(11,15),
Sigma = matrix(c(1,1/r,1/r,1), ncol = 2),
empirical = T)
capas[,,r] = z
}
capas[,1,85]
## [1] 9.929500 11.765780 11.682202 11.079675 10.422200 10.615726 11.240785
## [8] 10.779194 10.146963 10.876085 12.037225 10.143630 11.192545 11.195521
## [15] 10.825790 9.243551 11.446796 9.217345 10.924701 11.142292 11.859196
## [22] 11.379339 10.654452 10.062452 10.407207 10.456460 9.510922 11.060575
## [29] 12.384574 11.857318 9.558687 11.027170 10.361486 13.009447 10.536715
## [36] 11.134648 11.762018 12.079429 11.334655 11.537158 8.820859 9.726563
## [43] 13.699308 10.392943 13.347960 11.116443 10.826524 10.949695 11.167978
## [50] 12.070313
capas[,2,85]
## [1] 14.48606 13.67538 16.11193 15.04683 16.57760 13.37580 14.55400
## [8] 15.10468 15.03449 13.70344 14.59463 15.61050 15.42176 15.24519
## [15] 17.34281 14.59772 16.83544 14.04497 15.65233 16.04535 15.15141
## [22] 15.06630 13.67311 13.86030 15.15816 16.40232 14.92806 16.12617
## [29] 15.10059 15.31918 14.96830 14.83932 13.66117 14.20887 14.06237
## [36] 15.95870 15.10366 15.20313 14.24361 13.55827 16.60000 14.48370
## [43] 14.99518 12.67273 15.05247 15.38116 15.12892 15.19339 16.88428
## [50] 13.95427
plot(capas[,1,85], capas[,2,85])
cap = data.frame(cea75 = capas[,1,85], cea150 = capas[,2,85])
(cea75 = sort.int(rnorm(50, 11, 1), partial = 25))
## [1] 9.323085 9.097547 9.365965 9.788991 9.852177 9.403161 9.542438
## [8] 9.612030 9.722334 9.826791 9.786836 10.031096 10.273933 10.454872
## [15] 10.026384 10.449437 10.472958 10.528212 10.557556 10.607054 10.606804
## [22] 10.586087 10.567277 10.579762 10.743092 10.790699 10.814185 10.839288
## [29] 10.846451 10.961425 11.351966 11.214055 10.961730 11.527919 11.454708
## [36] 11.022761 11.120930 11.655514 11.432205 11.498524 11.354744 11.880011
## [43] 11.885862 12.150864 11.571159 11.882304 11.664067 13.093863 11.841346
## [50] 11.543646
(cea150 = sort.int(rnorm(50, 15, 1.2), partial = 25))
## [1] 12.13511 13.34864 13.44024 12.57692 12.26907 13.01241 13.67471
## [8] 13.74144 13.88663 13.87391 14.01893 13.69731 13.90118 14.06066
## [15] 14.08641 14.12283 14.21825 14.39100 14.41788 14.68973 14.79308
## [22] 14.78774 14.64356 14.61989 14.85132 15.00853 15.09480 16.19474
## [29] 16.03852 15.10489 15.80552 15.87487 16.23751 16.29508 15.87582
## [36] 15.67881 15.78143 16.01889 16.06085 16.04376 15.86776 15.79997
## [43] 16.15366 16.39186 16.30749 16.43007 16.91897 16.60619 16.42188
## [50] 17.62719
cap = data.frame(cea75, cea150)
cor(cap)
## cea75 cea150
## cea75 1.0000000 0.8948687
## cea150 0.8948687 1.0000000
(S = cov(cap))
## cea75 cea150
## cea75 0.7427641 0.9883052
## cea150 0.9883052 1.6421492
eigen(S)
## eigen() decomposition
## $values
## [1] 2.2782608 0.1066525
##
## $vectors
## [,1] [,2]
## [1,] 0.5412226 -0.8408794
## [2,] 0.8408794 0.5412226
1.853/(1.853+0.362)
## [1] 0.8365688
pca = princomp(cap)
pca$loadings
##
## Loadings:
## Comp.1 Comp.2
## cea75 0.541 0.841
## cea150 0.841 -0.541
##
## Comp.1 Comp.2
## SS loadings 1.0 1.0
## Proportion Var 0.5 0.5
## Cumulative Var 0.5 1.0
pca$scores
## Comp.1 Comp.2
## [1,] -3.16998054 0.327547845
## [2,] -2.27161884 -0.518889957
## [3,] -2.04931858 -0.342759700
## [4,] -2.54631162 0.480199797
## [5,] -2.77097863 0.699947291
## [6,] -2.38893574 -0.079933903
## [7,] -1.75664713 -0.321267281
## [8,] -1.66286968 -0.298864875
## [9,] -1.48107752 -0.284696879
## [10,] -1.43524319 -0.189974329
## [11,] -1.33492615 -0.302057351
## [12,] -1.47317088 0.077403875
## [13,] -1.17030926 0.171259341
## [14,] -0.93828046 0.237095469
## [15,] -1.14853275 -0.137149740
## [16,] -0.88894186 0.198875945
## [17,] -0.79597518 0.167010756
## [18,] -0.62080504 0.119974330
## [19,] -0.58231929 0.130099332
## [20,] -0.32693819 0.024591316
## [21,] -0.24016941 -0.031554381
## [22,] -0.25587348 -0.046083696
## [23,] -0.38729504 0.016134597
## [24,] -0.40043638 0.039440985
## [25,] -0.11743781 0.051529442
## [26,] 0.04052326 0.006475453
## [27,] 0.12577440 -0.020465186
## [28,] 1.06428071 -0.594671323
## [29,] 0.93679330 -0.504096469
## [30,] 0.21394952 0.097883695
## [31,] 1.01446933 0.047082315
## [32,] 0.99813801 -0.106415102
## [33,] 1.16651055 -0.514858821
## [34,] 1.52135687 -0.069922247
## [35,] 1.12918823 0.095428522
## [36,] 0.72974968 -0.161161766
## [37,] 0.86916869 -0.134151806
## [38,] 1.35817680 0.186847069
## [39,] 1.27259368 -0.023634315
## [40,] 1.29411891 0.041379368
## [41,] 1.06830819 0.015733064
## [42,] 1.29559278 0.494107337
## [43,] 1.59616374 0.307607111
## [44,] 1.93988964 0.401519741
## [45,] 1.55519465 -0.040278421
## [46,] 1.82667319 0.155010227
## [47,] 2.11966216 -0.293103247
## [48,] 2.63048960 1.078467237
## [49,] 1.79761385 0.125005239
## [50,] 2.65001290 -0.777665901
library(factoextra)
## Warning: package 'factoextra' was built under R version 3.6.1
## Loading required package: ggplot2
## Registered S3 methods overwritten by 'ggplot2':
## method from
## [.quosures rlang
## c.quosures rlang
## print.quosures rlang
## Welcome! Related Books: `Practical Guide To Cluster Analysis in R` at https://goo.gl/13EFCZ

# factoextra::fviz_pca(cap)
pca = prcomp(cap, scale = F)
fviz_pca(pca)

get_pca(pca)
## Principal Component Analysis Results for variables
## ===================================================
## Name Description
## 1 "$coord" "Coordinates for the variables"
## 2 "$cor" "Correlations between variables and dimensions"
## 3 "$cos2" "Cos2 for the variables"
## 4 "$contrib" "contributions of the variables"
#####
A = rnorm(100, 60, 5)
a = rnorm(100, 20, 2)
L = 100-(A+a)
plot(A,a)

library(rgl)
plot3d(A,a,L)
tex = data.frame(A,a,L)
princomp(tex)
## Call:
## princomp(x = tex)
##
## Standard deviations:
## Comp.1 Comp.2 Comp.3
## 7.583101e+00 2.411862e+00 6.977391e-08
##
## 3 variables and 100 observations.
eigen(cov(tex))
## eigen() decomposition
## $values
## [1] 5.808426e+01 5.875838e+00 -3.844020e-15
##
## $vectors
## [,1] [,2] [,3]
## [1,] 0.69529010 0.4280635 0.5773503
## [2,] 0.02306881 -0.8161706 0.5773503
## [3,] -0.71835890 0.3881071 0.5773503