実データに対する因子分析②
南山大学 浦上昌則 先生のRの練習帳で学習した記録です。この場を借りて御礼申し上げます。
library(psych)
library(MASS)
#3D表示
options(rgl.printRglwidget = TRUE)
#WebGLに変換してRStudio, VSCode, ウェブブラウザなど出力
library(rgl)
library(ggplot2)##
## Attaching package: 'ggplot2'## The following objects are masked from 'package:psych':
##
## %+%, alphalibrary(ggExtra)- 第 1 ステップは，固有値を求めて，抽出因子数の目処を立てること。
1 データ
x <- read.csv('urakami_data1.csv',header=TRUE)
#head(x)
#x[5:7]#列選択
#cor(x[5:7])
#cor(x[5:7], use="complete.obs")#リストワイズ削除
#cor(x[5:7], use="pairwise.complete.obs")#2変数の組合せごとに欠損値を含むサンプルを取り除いて計算
label_b <- c("b1","b2","b3","b4","b5","b6","b7","b8","b9","b10","b11",
"b12","b13","b14","b15","b16","b17","b18","b19","b20")
xb <- x[label_b]#必要な列を抽出して，新しい名前でまとめておく
#xb
#cor(xb, use="complete.obs")
#round(cor(xb, use="complete.obs"),3)2 因子分析(1)
因子数を探索
library("psych")
VSS.scree(xb)#スクリープロット
xb.cor <- cor(xb, use="complete.obs") #NAを含む組み合わせを取り除いたデータで相関係数
eigen(xb.cor)$value#eigen$valueを使う方法もある## [1] 6.9395390 2.3645976 1.9234457 1.2390728 1.1241792 1.0556399 0.6711797
## [8] 0.5972908 0.4926267 0.4485218 0.4391389 0.4099862 0.3649860 0.3338334
## [15] 0.3307298 0.2916840 0.2782851 0.2569604 0.2267127 0.2115902plot(eigen(xb.cor)$value,type="b",main="Scree Plot", xlab="NUM", ylab="VALUE")
fa.parallel(xb)#平行分析
## Parallel analysis suggests that the number of factors = 6 and the number of components = 3平行分析
分析するデータと同じ変数の数，同じサンプル数の乱数行列を作り，そこから相関行列を導き固有値を求めるという作業を行っています。
ランダムな状況から生じた固有値の推移と分析データの推移を比較し，ランダムなものよりも大きい固有値の個数だけ因子を抽出しようという考え方です。
出力される図の中で赤の波線になっている部分が，乱数から作られた固有値の推移です。
3 因子分析
library(GPArotation)##
## Attaching package: 'GPArotation'## The following objects are masked from 'package:psych':
##
## equamax, variminf1 <- fa(xb, nfactors=3, fm="pa", rotate="promax") #fm=は因子抽出方法、paは主因子法
#"minres"ミンレス法(最小残差法)，"pa"主因子法，"ml"て再尤法，"gls"一般化最小二乗法
print(f1, sort=TRUE, digit=3)## Factor Analysis using method = pa
## Call: fa(r = xb, nfactors = 3, rotate = "promax", fm = "pa")
## Standardized loadings (pattern matrix) based upon correlation matrix
## item PA1 PA2 PA3 h2 u2 com
## b5 5 0.796 0.215 -0.118 0.734 0.266 1.19
## b14 14 0.773 0.135 -0.069 0.653 0.347 1.08
## b20 20 -0.729 -0.234 0.321 0.555 0.445 1.60
## b13 13 0.664 -0.072 -0.184 0.319 0.681 1.18
## b1 1 0.586 -0.148 -0.120 0.240 0.760 1.22
## b2 2 0.554 -0.219 0.220 0.392 0.608 1.64
## b16 16 0.464 0.044 0.223 0.398 0.602 1.46
## b10 10 0.447 -0.119 0.156 0.251 0.749 1.39
## b4 4 0.420 0.081 0.403 0.569 0.431 2.07
## b11 11 -0.207 0.883 -0.003 0.658 0.342 1.11
## b6 6 -0.051 0.820 -0.006 0.636 0.364 1.01
## b3 3 -0.166 0.697 0.024 0.417 0.583 1.12
## b18 18 0.052 0.570 -0.011 0.350 0.650 1.02
## b8 8 0.241 0.437 -0.053 0.320 0.680 1.59
## b15 15 0.198 0.423 0.374 0.604 0.396 2.41
## b19 19 -0.313 -0.028 0.811 0.489 0.511 1.29
## b9 9 -0.176 -0.025 0.757 0.460 0.540 1.11
## b12 12 0.366 0.032 0.499 0.593 0.407 1.84
## b17 17 0.356 0.115 0.469 0.601 0.399 2.00
## b7 7 0.421 0.008 0.421 0.542 0.458 2.00
##
## PA1 PA2 PA3
## SS loadings 4.291 2.948 2.542
## Proportion Var 0.215 0.147 0.127
## Cumulative Var 0.215 0.362 0.489
## Proportion Explained 0.439 0.301 0.260
## Cumulative Proportion 0.439 0.740 1.000
##
## With factor correlations of
## PA1 PA2 PA3
## PA1 1.000 0.446 0.515
## PA2 0.446 1.000 0.301
## PA3 0.515 0.301 1.000
##
## Mean item complexity = 1.5
## Test of the hypothesis that 3 factors are sufficient.
##
## The degrees of freedom for the null model are 190 and the objective function was 10.449 with Chi Square of 3025.077
## The degrees of freedom for the model are 133 and the objective function was 2.004
##
## The root mean square of the residuals (RMSR) is 0.06
## The df corrected root mean square of the residuals is 0.072
##
## The harmonic number of observations is 298 with the empirical chi square 407.646 with prob < 1.2e-29
## The total number of observations was 298 with Likelihood Chi Square = 576.008 with prob < 1.99e-56
##
## Tucker Lewis Index of factoring reliability = 0.7751
## RMSEA index = 0.1057 and the 90 % confidence intervals are 0.0971 0.1149
## BIC = -181.706
## Fit based upon off diagonal values = 0.969
## Measures of factor score adequacy
## PA1 PA2 PA3
## Correlation of (regression) scores with factors 0.954 0.936 0.920
## Multiple R square of scores with factors 0.910 0.875 0.847
## Minimum correlation of possible factor scores 0.819 0.750 0.694biplot(f1)
fa.diagram( f1 )
3.1 ３次元にして図示
d <- data.frame(X=xb[,1],Y=xb[,2],Z=xb[,3])
attach(d)
str(d)## 'data.frame': 298 obs. of 3 variables:
## $ X: int 4 4 NA 2 3 4 2 4 4 4 ...
## $ Y: int 4 3 4 4 2 4 2 3 4 3 ...
## $ Z: int 2 2 1 1 3 1 1 4 3 2 ...xm <- mean(X,na.rm=T)
ym <- mean(Y)
zm <- mean(Z)
fit <- fa(r=d, nfactors=2 ,rotate="", fm="ml", scores=T)## Specified rotation not found, rotate='none' usedprint(fit, digits=5,sort=T)# 因子得点を出すためにはscores=T## Factor Analysis using method = ml
## Call: fa(r = d, nfactors = 2, rotate = "", scores = T, fm = "ml")
## Standardized loadings (pattern matrix) based upon correlation matrix
## item ML1 ML2 h2 u2 com
## X 1 0.65228 0.26416 0.495256 0.50474 1.3194
## Y 2 0.58418 -0.33152 0.451174 0.54883 1.5836
## Z 3 0.09971 0.11266 0.022634 0.97737 1.9709
##
## ML1 ML2
## SS loadings 0.77669 0.19238
## Proportion Var 0.25890 0.06413
## Cumulative Var 0.25890 0.32302
## Proportion Explained 0.80148 0.19852
## Cumulative Proportion 0.80148 1.00000
##
## Mean item complexity = 1.6
## Test of the hypothesis that 2 factors are sufficient.
##
## The degrees of freedom for the null model are 3 and the objective function was 0.09915 with Chi Square of 29.26469
## The degrees of freedom for the model are -2 and the objective function was 0
##
## The root mean square of the residuals (RMSR) is 0
## The df corrected root mean square of the residuals is NA
##
## The harmonic number of observations is 297 with the empirical chi square 0 with prob < NA
## The total number of observations was 298 with Likelihood Chi Square = 0 with prob < NA
##
## Tucker Lewis Index of factoring reliability = 1.1148
## Fit based upon off diagonal values = 1
## Measures of factor score adequacy
## ML1 ML2
## Correlation of (regression) scores with factors 0.77198 0.50998
## Multiple R square of scores with factors 0.59595 0.26008
## Minimum correlation of possible factor scores 0.19190 -0.47985biplot(fit)
fa.diagram( fit )
plot3d(X, Y, Z, type = "s", col = "gray", size = 1)
v <- fit$ loadings
#v[1,1];v[2,1];v[3,1]
#v[1,2];v[2,2];v[3,2]
#str(v)
arrow3d(p0 = c(xm,ym,zm), p1 = c(xm,ym,zm)+c(v[1,1],v[2,1],v[3,1]), col = "red")
arrow3d(p0 = c(xm,ym,zm), p1 = c(xm,ym,zm)+c(v[1,2],v[2,2],v[3,2]), col = "blue")3.2 ２次元
plot(jitter(X),jitter(Y))
arrows(x0 = xm, y0 = ym, x1 = xm+v[1,1], y1 =ym+v[2,1], col = "red")
arrows(x0 = xm, y0 = ym, x1 = xm+v[1,2], y1 =ym+v[2,2], col = "blue")
4 因子分析(2)
4.1 プロマックス回転
library("psych")
f1 <- fa(xb, nfactors=3, fm="ml", rotate="promax") #fm=は因子抽出方法、mlは最尤法
print(f1, sort=TRUE, digit=3)## Factor Analysis using method = ml
## Call: fa(r = xb, nfactors = 3, rotate = "promax", fm = "ml")
## Standardized loadings (pattern matrix) based upon correlation matrix
## item ML1 ML3 ML2 h2 u2 com
## b20 20 -0.800 0.228 -0.162 0.604 0.396 1.25
## b14 14 0.785 0.026 0.070 0.691 0.309 1.02
## b5 5 0.763 0.017 0.174 0.737 0.263 1.10
## b13 13 0.633 -0.092 -0.080 0.314 0.686 1.08
## b1 1 0.552 -0.060 -0.145 0.233 0.767 1.16
## b16 16 0.364 0.310 0.042 0.374 0.626 1.98
## b10 10 0.312 0.258 -0.102 0.220 0.780 2.17
## b12 12 0.165 0.690 0.027 0.644 0.356 1.12
## b17 17 0.165 0.668 0.096 0.654 0.346 1.17
## b19 19 -0.308 0.657 -0.033 0.302 0.698 1.42
## b9 9 -0.199 0.632 -0.013 0.299 0.701 1.20
## b7 7 0.249 0.597 -0.014 0.573 0.427 1.34
## b4 4 0.257 0.554 0.076 0.576 0.424 1.45
## b15 15 0.060 0.546 0.393 0.642 0.358 1.84
## b2 2 0.362 0.373 -0.182 0.356 0.644 2.44
## b11 11 -0.174 -0.010 0.903 0.714 0.286 1.07
## b6 6 -0.011 0.035 0.790 0.635 0.365 1.00
## b3 3 -0.158 0.005 0.752 0.495 0.505 1.09
## b18 18 0.133 -0.011 0.494 0.311 0.689 1.15
## b8 8 0.351 -0.064 0.351 0.313 0.687 2.07
##
## ML1 ML3 ML2
## SS loadings 3.610 3.305 2.773
## Proportion Var 0.181 0.165 0.139
## Cumulative Var 0.181 0.346 0.484
## Proportion Explained 0.373 0.341 0.286
## Cumulative Proportion 0.373 0.714 1.000
##
## With factor correlations of
## ML1 ML3 ML2
## ML1 1.000 0.545 0.410
## ML3 0.545 1.000 0.304
## ML2 0.410 0.304 1.000
##
## Mean item complexity = 1.4
## Test of the hypothesis that 3 factors are sufficient.
##
## The degrees of freedom for the null model are 190 and the objective function was 10.449 with Chi Square of 3025.077
## The degrees of freedom for the model are 133 and the objective function was 1.935
##
## The root mean square of the residuals (RMSR) is 0.063
## The df corrected root mean square of the residuals is 0.075
##
## The harmonic number of observations is 298 with the empirical chi square 442.51 with prob < 6.71e-35
## The total number of observations was 298 with Likelihood Chi Square = 556.348 with prob < 3.84e-53
##
## Tucker Lewis Index of factoring reliability = 0.7851
## RMSEA index = 0.1033 and the 90 % confidence intervals are 0.0947 0.1125
## BIC = -201.366
## Fit based upon off diagonal values = 0.966
## Measures of factor score adequacy
## ML1 ML3 ML2
## Correlation of (regression) scores with factors 0.948 0.937 0.937
## Multiple R square of scores with factors 0.899 0.878 0.878
## Minimum correlation of possible factor scores 0.799 0.757 0.755print( f1$loadings,sort=TRUE, digits = 2, cutoff = 0.3 )##
## Loadings:
## ML1 ML3 ML2
## b1 0.55
## b5 0.76
## b13 0.63
## b14 0.79
## b20 -0.80
## b4 0.55
## b7 0.60
## b9 0.63
## b12 0.69
## b15 0.55 0.39
## b17 0.67
## b19 -0.31 0.66
## b3 0.75
## b6 0.79
## b11 0.90
## b2 0.36 0.37
## b8 0.35 0.35
## b10 0.31
## b16 0.36 0.31
## b18 0.49
##
## ML1 ML3 ML2
## SS loadings 3.42 3.09 2.68
## Proportion Var 0.17 0.15 0.13
## Cumulative Var 0.17 0.33 0.46#comは複雑性biplot(f1)
fa.diagram( f1 )
回転の種類
4.2 独立クラスター回転
各項目に1因子しか負荷しない状態を目指す
f1 <- fa(xb, nfactors=3, fm="ml", rotate="cluster") #fm=は因子抽出方法、mlは最尤法
print(f1, sort=TRUE, cutoff=0.3)## Factor Analysis using method = ml
## Call: fa(r = xb, nfactors = 3, rotate = "cluster", fm = "ml")
## Standardized loadings (pattern matrix) based upon correlation matrix
## item ML3 ML1 ML2 h2 u2 com
## b12 12 0.73 0.14 -0.02 0.64 0.36 1.1
## b17 17 0.71 0.14 0.05 0.65 0.35 1.1
## b19 19 0.66 -0.32 -0.09 0.30 0.70 1.5
## b9 9 0.65 -0.21 -0.06 0.30 0.70 1.2
## b7 7 0.64 0.22 -0.06 0.57 0.43 1.3
## b4 4 0.59 0.23 0.04 0.58 0.42 1.3
## b15 15 0.58 0.04 0.35 0.64 0.36 1.7
## b2 2 0.41 0.33 -0.20 0.36 0.64 2.4
## b16 16 0.35 0.34 0.02 0.37 0.63 2.0
## b20 20 0.19 -0.78 -0.19 0.60 0.40 1.2
## b14 14 0.08 0.75 0.08 0.69 0.31 1.0
## b5 5 0.07 0.73 0.18 0.74 0.26 1.1
## b13 13 -0.06 0.61 -0.06 0.31 0.69 1.0
## b1 1 -0.03 0.53 -0.13 0.23 0.77 1.1
## b10 10 0.29 0.29 -0.12 0.22 0.78 2.3
## b11 11 -0.01 -0.16 0.89 0.71 0.29 1.1
## b6 6 0.04 -0.01 0.78 0.64 0.36 1.0
## b3 3 0.00 -0.15 0.74 0.49 0.51 1.1
## b18 18 0.00 0.13 0.49 0.31 0.69 1.1
## b8 8 -0.04 0.34 0.36 0.31 0.69 2.0
##
## ML3 ML1 ML2
## SS loadings 3.62 3.39 2.68
## Proportion Var 0.18 0.17 0.13
## Cumulative Var 0.18 0.35 0.48
## Proportion Explained 0.37 0.35 0.28
## Cumulative Proportion 0.37 0.72 1.00
##
## With factor correlations of
## ML3 ML1 ML2
## ML3 1.00 0.54 0.36
## ML1 0.54 1.00 0.38
## ML2 0.36 0.38 1.00
##
## Mean item complexity = 1.4
## Test of the hypothesis that 3 factors are sufficient.
##
## The degrees of freedom for the null model are 190 and the objective function was 10.45 0.3 with Chi Square of 3025.08
## The degrees of freedom for the model are 133 and the objective function was 1.94
## 0.3
## The root mean square of the residuals (RMSR) is 0.06
## The df corrected root mean square of the residuals is 0.07
## 0.3
## The harmonic number of observations is 298 with the empirical chi square 442.51 with prob < 6.7e-35
## 0.3The total number of observations was 298 with Likelihood Chi Square = 556.35 with prob < 3.8e-53
## 0.3
## Tucker Lewis Index of factoring reliability = 0.785
## RMSEA index = 0.103 and the 90 % confidence intervals are 0.095 0.113 0.3
## BIC = -201.37
## Fit based upon off diagonal values = 0.97
## Measures of factor score adequacy
## ML3 ML1 ML2
## Correlation of (regression) scores with factors 0.94 0.94 0.94
## Multiple R square of scores with factors 0.89 0.89 0.88
## Minimum correlation of possible factor scores 0.78 0.78 0.75print( f1$loadings,sort=TRUE, digits = 2, cutoff = 0.3 )##
## Loadings:
## ML3 ML1 ML2
## b4 0.59
## b7 0.64
## b9 0.65
## b12 0.73
## b15 0.58 0.35
## b17 0.71
## b19 0.66 -0.32
## b1 0.53
## b5 0.73
## b13 0.61
## b14 0.75
## b20 -0.78
## b3 0.74
## b6 0.78
## b11 0.89
## b2 0.41 0.33
## b8 0.34 0.36
## b10
## b16 0.35 0.34
## b18 0.49
##
## ML3 ML1 ML2
## SS loadings 3.41 3.13 2.63
## Proportion Var 0.17 0.16 0.13
## Cumulative Var 0.17 0.33 0.46biplot(f1)
fa.diagram( f1 )
label_f1 <- c("b12","b17","b7","b15","b4","b9")
label_f2 <- c("b20","b14","b5","b13","b1")
label_f3 <- c("b11","b3","b6","b18")
dat_f1 <- xb[label_f1]
dat_f2 <- xb[label_f2]
dat_f3 <- xb[label_f3]library(psych)
psych::alpha(dat_f1)##
## Reliability analysis
## Call: psych::alpha(x = dat_f1)
##
## raw_alpha std.alpha G6(smc) average_r S/N ase mean sd median_r
## 0.86 0.86 0.85 0.51 6.3 0.013 3.1 0.68 0.59
##
## 95% confidence boundaries
## lower alpha upper
## Feldt 0.83 0.86 0.88
## Duhachek 0.83 0.86 0.88
##
## Reliability if an item is dropped:
## raw_alpha std.alpha G6(smc) average_r S/N alpha se var.r med.r
## b12 0.82 0.83 0.81 0.49 4.7 0.016 0.0237 0.54
## b17 0.82 0.82 0.80 0.48 4.7 0.017 0.0213 0.54
## b7 0.83 0.83 0.82 0.49 4.9 0.016 0.0273 0.57
## b15 0.83 0.84 0.82 0.51 5.1 0.016 0.0225 0.60
## b4 0.83 0.83 0.82 0.49 4.9 0.016 0.0290 0.57
## b9 0.89 0.89 0.87 0.61 7.8 0.010 0.0021 0.61
##
## Item statistics
## n raw.r std.r r.cor r.drop mean sd
## b12 298 0.82 0.82 0.79 0.73 3.3 0.84
## b17 298 0.84 0.83 0.81 0.75 3.0 0.93
## b7 298 0.80 0.81 0.77 0.71 3.4 0.83
## b15 298 0.79 0.78 0.74 0.67 2.6 0.98
## b4 298 0.80 0.81 0.76 0.71 3.4 0.84
## b9 298 0.56 0.56 0.41 0.38 2.8 0.93
##
## Non missing response frequency for each item
## 1 2 3 4 miss
## b12 0.05 0.11 0.37 0.47 0
## b17 0.08 0.18 0.39 0.35 0
## b7 0.05 0.07 0.31 0.57 0
## b15 0.16 0.29 0.36 0.20 0
## b4 0.05 0.07 0.31 0.56 0
## b9 0.11 0.27 0.40 0.23 0psych::alpha(dat_f2)#エラー回避## Warning in psych::alpha(dat_f2): Some items were negatively correlated with the total scale and probably
## should be reversed.
## To do this, run the function again with the 'check.keys=TRUE' option## Some items ( b20 ) were negatively correlated with the total scale and
## probably should be reversed.
## To do this, run the function again with the 'check.keys=TRUE' option##
## Reliability analysis
## Call: psych::alpha(x = dat_f2)
##
## raw_alpha std.alpha G6(smc) average_r S/N ase mean sd median_r
## 0.18 0.32 0.6 0.087 0.48 0.067 2.6 0.49 0.38
##
## 95% confidence boundaries
## lower alpha upper
## Feldt 0.03 0.18 0.32
## Duhachek 0.05 0.18 0.31
##
## Reliability if an item is dropped:
## raw_alpha std.alpha G6(smc) average_r S/N alpha se var.r med.r
## b20 0.78 0.782 0.76 0.4729 3.589 0.021 0.015 0.446
## b14 -0.19 0.012 0.38 0.0030 0.012 0.110 0.246 0.072
## b5 -0.26 -0.030 0.37 -0.0072 -0.029 0.117 0.246 0.052
## b13 -0.26 -0.058 0.43 -0.0139 -0.055 0.112 0.338 0.052
## b1 -0.34 -0.078 0.44 -0.0184 -0.072 0.116 0.380 0.015
##
## Item statistics
## n raw.r std.r r.cor r.drop mean sd
## b20 298 -0.31 -0.37 -0.73 -0.64 2.3 1.17
## b14 298 0.69 0.71 0.72 0.37 2.7 0.97
## b5 298 0.72 0.74 0.74 0.40 2.9 1.04
## b13 298 0.73 0.75 0.65 0.46 2.5 0.90
## b1 297 0.76 0.76 0.63 0.49 2.6 0.97
##
## Non missing response frequency for each item
## 1 2 3 4 miss
## b20 0.36 0.18 0.24 0.21 0
## b14 0.14 0.26 0.38 0.21 0
## b5 0.14 0.19 0.32 0.35 0
## b13 0.15 0.37 0.35 0.13 0
## b1 0.14 0.29 0.36 0.20 0psych::alpha(dat_f3)##
## Reliability analysis
## Call: psych::alpha(x = dat_f3)
##
## raw_alpha std.alpha G6(smc) average_r S/N ase mean sd median_r
## 0.8 0.8 0.77 0.49 3.9 0.018 2 0.8 0.47
##
## 95% confidence boundaries
## lower alpha upper
## Feldt 0.76 0.8 0.84
## Duhachek 0.77 0.8 0.84
##
## Reliability if an item is dropped:
## raw_alpha std.alpha G6(smc) average_r S/N alpha se var.r med.r
## b11 0.68 0.68 0.60 0.41 2.1 0.030 0.0108 0.43
## b3 0.75 0.74 0.68 0.49 2.9 0.024 0.0203 0.43
## b6 0.72 0.72 0.67 0.46 2.5 0.026 0.0368 0.39
## b18 0.82 0.83 0.77 0.61 4.8 0.018 0.0077 0.66
##
## Item statistics
## n raw.r std.r r.cor r.drop mean sd
## b11 298 0.87 0.86 0.84 0.75 2.3 1.0
## b3 298 0.80 0.79 0.71 0.62 2.0 1.1
## b6 298 0.84 0.82 0.75 0.67 2.2 1.1
## b18 298 0.63 0.67 0.48 0.44 1.6 0.8
##
## Non missing response frequency for each item
## 1 2 3 4 miss
## b11 0.29 0.31 0.26 0.14 0
## b3 0.42 0.29 0.16 0.13 0
## b6 0.40 0.19 0.23 0.17 0
## b18 0.52 0.35 0.09 0.04 0cor(dat_f1, use="complete.obs")## b12 b17 b7 b15 b4 b9
## b12 1.0000000 0.6492228 0.6265161 0.6131552 0.6130023 0.3100640
## b17 0.6492228 1.0000000 0.6396938 0.6799454 0.5860129 0.3037255
## b7 0.6265161 0.6396938 1.0000000 0.5210908 0.6064429 0.3561818
## b15 0.6131552 0.6799454 0.5210908 1.0000000 0.5611867 0.2494030
## b4 0.6130023 0.5860129 0.6064429 0.5611867 1.0000000 0.3664278
## b9 0.3100640 0.3037255 0.3561818 0.2494030 0.3664278 1.0000000cor.plot(cor(dat_f1, use="complete.obs"), numbers=TRUE)