Nothing
# file OrdinalLogisticBiplot/R/OrdinalLogisticBiplot.R
# copyright (C) 2012-2013 J.C. Hernandez and J.L. Vicente-Villardon
#
# This program is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 2 or 3 of the License
# (at your option).
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# A copy of the GNU General Public License is available at
# http://www.r-project.org/Licenses/
#
OrdinalLogisticBiplot <- function(datanom,sFormula=NULL,numFactors=2,method="EM",rotation="varimax",
metfsco="EAP",nnodos = 10, tol = 1e-04, maxiter = 100, penalization = 0.1,cte=TRUE,
show=FALSE,ItemCurves = FALSE,initial=1,alfa=1){
#We have to check if datanom is a data frame or a matrix to tackle with sFormula parameter
if(!(is.null(sFormula))){
if(is.data.frame(datanom)){
datanom = model.frame(formula=sFormula,data=datanom)
}else if(is.matrix(datanom)){
datanom = model.frame(formula=sFormula,data=as.data.frame(datanom))
datanom = as.matrix(datanom)
}else{
print("It is not posible to use the formula passed as parameter. Data are not a data frame nor matrix")
}
}
if(ncol(datanom) <= numFactors){
stop("It is not posible to reduce dimension because number of Factors is lower than variables in our data set")
}
dataSet = CheckDataSet(datanom)
#We have now in dataSet structure all our data, the names for the row and colums and de levels if exists
datanom = dataSet$datanom
nRowsdata <- dim(datanom)[1]
nColsdata <- dim(datanom)[2]
numVar <- ncol(datanom)
datanomcats = apply(datanom[,1:numVar], 2, function(x) nlevels(as.factor(x)))
x <- matrix(0,nRowsdata,numFactors)
if (method == "EM"){
olb = OrdinalLogBiplotEM(datanom,dim = numFactors, nnodos = nnodos, tol = tol, maxiter = maxiter,
penalization = penalization,initial=initial,show=show,alfa=alfa)
print(olb$ColumnParameters$fit)
coefs = olb$ColumnParameters$coefficients
thresholds = olb$ColumnParameters$thresholds
x = as.matrix(olb$RowCoordinates)
Fitting = olb$ColumnParameters$fit
LogLik = olb$LogLikelihood
FactorLoadings = olb$loadings
h2 = olb$r2
rotation=" "
metfsco=" "
}else if (method == "MIRT"){
rows = EstimationRowsMIRT(datanom,numFactors = numFactors,metfsco=metfsco,rotation=rotation,maxiter=maxiter) #Se estima una vez y luego se elige el plano que queremos estudiar en GetOrdinalBiplotObjectMIRT
fittedMirtModels = CalculateFittingIndicatorsMirt(rows)
x=as.matrix(rows$estimRows)
coefs = as.matrix(rows$sepCoefMirt$xCoeffic)
thresholds = rows$sepCoefMirt$indCoeffic
fittedInd = matrix(0,nColsdata,8)
dimnames(fittedInd)[[1]]= dimnames(datanom)[[2]][1:dim(datanom)[2]]
dimnames(fittedInd)[[2]]= c("logLik","Deviance","df","p-value","PCC","CoxSnell","Macfaden","Nagelkerke")
sumLogLik = 0
for(i in 1:nColsdata){
fittedInd[i,]=c(fittedMirtModels[,i]$logLik,fittedMirtModels[,i]$Deviance,
fittedMirtModels[,i]$df,fittedMirtModels[,i]$pval,
fittedMirtModels[,i]$PercentClasif,fittedMirtModels[,i]$CoxSnell,
fittedMirtModels[,i]$MacFaden,fittedMirtModels[,i]$Nagelkerke)
sumLogLik = sumLogLik + fittedMirtModels[,i]$logLik
}
Fitting = fittedInd
LogLik = sumLogLik
if(numFactors > 1){
FactorLoadings = rows$summ$rotF
h2 = rows$summ$h2
}else{
FactorLoadings = rows$summ$F
h2 = ""
}
rotation=rotation
metfsco=metfsco
nnodos=" "
olb = rows
}
if((numFactors > 1)|(method == "EM")){
FactorLoadingsComm = cbind(FactorLoadings,h2)
dimnames(FactorLoadingsComm)[[2]] = c(paste(c(rep("F_",numFactors)),c(1:numFactors),sep=""),"Communalities")
}else{
FactorLoadingsComm = FactorLoadings
dimnames(FactorLoadingsComm)[[2]] = "F_1"
}
dimnames(FactorLoadingsComm)[[1]] = dataSet$ColumNames
dimnames(x)[[1]]=dimnames(datanom)[[1]]
dimnames(x)[[2]]=c(1:numFactors)
dimnames(coefs)[[1]] = dataSet$ColumNames
dimnames(coefs)[[2]] = paste(c(rep("Dim",numFactors)),c(1:numFactors))
dimnames(thresholds)[[1]] = dataSet$ColumNames
dimnames(thresholds)[[2]] = c(1:ncol(thresholds))
ordinal.logistic.biplot<-list()
ordinal.logistic.biplot$dataSet = dataSet
ordinal.logistic.biplot$RowCoords = x
ordinal.logistic.biplot$Ncats = datanomcats
ordinal.logistic.biplot$estimObject = olb
ordinal.logistic.biplot$Fitting = Fitting
ordinal.logistic.biplot$coefs = coefs
ordinal.logistic.biplot$thresholds = thresholds
ordinal.logistic.biplot$NumFactors = numFactors
ordinal.logistic.biplot$Coordinates = method
ordinal.logistic.biplot$Rotation = rotation
ordinal.logistic.biplot$Methodfscores = metfsco
ordinal.logistic.biplot$NumNodos = nnodos
ordinal.logistic.biplot$tol = tol
ordinal.logistic.biplot$maxiter = maxiter
ordinal.logistic.biplot$penalization = penalization
ordinal.logistic.biplot$cte = cte
ordinal.logistic.biplot$show = show
ordinal.logistic.biplot$ItemCurves = ItemCurves
ordinal.logistic.biplot$LogLik = LogLik
ordinal.logistic.biplot$FactorLoadingsComm = FactorLoadingsComm
class(ordinal.logistic.biplot)='ordinal.logistic.biplot'
return(ordinal.logistic.biplot)
}
#This function shows a summary of the principal results from the estimation for rows and variables.
#----------------------Parameters
#x: object with the information needed about all the variables and individuals
summary.ordinal.logistic.biplot <- function(object,data = FALSE,rowCoords = FALSE,coefs = FALSE,loadCommun = FALSE,...) {
x=object
if(x$Coordinates == "EM"){
cat(paste(" Ordinal Logistic Biplot Estimation ", "with Ridge Penalization :", x$penalization, ", EM algorithm and logit link"), "\n")
cat("\n Percentage of correct classifications,Pseudo R-squared measures and other indicators: \n")
print(x$Fitting)
}else{
cat(paste(" Ordinal Logistic Biplot Estimation ", "using MIRT method"), "\n")
cat(paste("Rotation", x$Rotation), "\n")
cat(paste("Method of fscores calculation:", x$Methodfscores), "\n")
cat("\n Percentage of correct classifications,Pseudo R-squared measures and other indicators: \n")
print(x$Fitting)
}
cat(paste("\n Number of factors for the reduced solution:",x$NumFactors,"\n",sep=""))
cat("\n Number of categories of the variables:\n")
print(x$Ncats)
if(data){
cat("\n n: ", nrow(x$dataSet$datanom), "\n")
}
if(rowCoords){
cat("\nCoordinates for the individuals: ","\n")
print(x$RowCoords)
}
if(coefs){
cat("\n Coefficients:\n")
print(x$coefs)
cat("\n Thresholds:\n")
print(x$thresholds)
}
if(loadCommun){
cat("\n Factor Loadings and Communalities:\n")
print(x$FactorLoadingsComm)
}
}
plot.ordinal.logistic.biplot <- function(x,planex=1,planey=2,AtLeastR2 = 0.01,
xlimi=-1.5,xlimu=1.5,ylimi=-1.5,ylimu=1.5,margin = 0,
ShowAxis = TRUE, PlotVars = TRUE, PlotInd = TRUE, LabelVar = TRUE,
LabelInd = TRUE,CexInd = NULL, CexVar = NULL, ColorInd = NULL, ColorVar = NULL,
PchInd = NULL, PchVar = NULL,showIIC=FALSE,iicxi=-1.5,iicxu=1.5,
legendPlot = FALSE,PlotClus = FALSE,Clusters=NULL,chulls = TRUE,centers = TRUE,
colorCluster = NULL,ConfidentLevel=NULL,addToExistingPlot=FALSE,...) {
olbo = x
if(olbo$NumFactors == 1){
stop("There is only one factor and it is not posible to represent the biplot on two dimensions")
}
if(planex == planey){
stop("The plane of the biplot is not correct. It should be selected different factors.")
}
if((planex > olbo$NumFactors) | (planey > olbo$NumFactors)){
stop(paste("With ",olbo$NumFactors, " factors it is not posible to analyze the plane ",planex,"-",planey,". Please,
select another plane.",sep=""))
}
n = nrow(olbo$dataSet$datanom)
p = ncol(olbo$dataSet$datanom)
RowNames = olbo$dataSet$RowNames
VarNames = olbo$dataSet$ColumNames
DimNames = "Dim 1"
for (i in 2:olbo$NumFactors)
DimNames = c(DimNames, paste("Dim", i))
# Determining sizes and colors of the points
if (is.null(CexInd)){
CexInd = rep(0.5, n)
}else{
if (length(CexInd) == 1){
CexInd = rep(CexInd, n)
}else if(length(CexInd) < n){
CexInd = rep(0.5, n)
}else{
CexInd = CexInd[1:n]
}
}
if (is.null(ColorInd)){
ColorInd = rep("black", n)
}else{
if (length(ColorInd) == 1){
ColorInd = rep(ColorInd, n)
}else if(length(ColorInd) < n){
ColorInd = rep("black", n)
}else{
ColorInd = ColorInd[1:n]
}
}
if (is.null(PchInd)){
PchInd = rep(1, n)
}else{
if (length(PchInd) == 1){
PchInd = rep(PchInd, n)
}else if(length(PchInd) < n){
PchInd = rep(1, n)
}else{
PchInd = PchInd[1:n]
}
}
if (is.null(CexVar)){
CexVar = rep(0.8, p)
}else{
if (length(CexVar) == 1){
CexVar = rep(CexVar, p)
}else if(length(CexVar) < p){
print("It has been specified lower cex values for the variables than variables")
CexVar = rep(0.8, p)
}else{
CexVar = CexVar[1:p]
}
}
if (is.null(PchVar)){
PchVar = c(0:(p-1))
}else{
if (length(PchVar) == 1){
PchVar = c(0:(p-1))
}else if(length(PchVar) < p){
print("It has been specified lower pch values for the variables than variables")
PchVar = c(0:(p-1))
}else{
PchVar = PchVar[1:p]
}
}
if (is.null(ColorVar)){
ColorVar = colors()[20 + 2*c(1:p)]
}else{
if (length(ColorVar) == 1){
ColorVar = colors()[20 + 2*c(1:p)]
}else if(length(ColorVar) < p){
print("It has been specified lower color values for the variables than variables")
ColorVar = colors()[20 + 2*c(1:p)]
}else{
ColorVar = ColorVar[1:p]
}
}
if (ShowAxis) {
xaxt = "s"
yaxt = "s"
} else {
xaxt = "n"
yaxt = "n"
}
if ((margin < 0) | (margin > 0.3))
margin = 0
xmin= xlimi - (xlimu - xlimi) * margin
xmax= xlimu + (xlimu - xlimi) * margin
ymin= ylimi - (ylimu - ylimi) * margin
ymax= ylimu + (ylimu - ylimi) * margin
if(PlotInd == TRUE){
if(addToExistingPlot == TRUE){
points(olbo$RowCoords[,planex], olbo$RowCoords[,planey], pch = PchInd, col=ColorInd, cex = CexInd)
}else{
dev.new()
plot(olbo$RowCoords[,planex], olbo$RowCoords[,planey], cex = CexInd, col=ColorInd, pch = PchInd, asp=1, xaxt = xaxt, yaxt = yaxt ,xlim=c(xmin,xmax),ylim=c(ymin,ymax),
main="Ordinal Logistic Biplot", xlab=paste("Axis ",planex,sep=""), ylab=paste("Axis ",planey,sep=""))
}
if(LabelInd == TRUE){
text(olbo$RowCoords[,planex], olbo$RowCoords[,planey],row.names(RowNames), cex = CexInd,col=ColorInd,pos=1,offset=0.1)
}
}else{
if(!addToExistingPlot){
plot(olbo$RowCoords[,planex], olbo$RowCoords[,planey], cex = 0,asp=1, xaxt = xaxt, yaxt = yaxt ,xlim=c(xmin,xmax),ylim=c(ymin,ymax),
main="Ordinal Logistic Biplot", xlab=paste("Axis ",planex,sep=""), ylab=paste("Axis ",planey,sep=""))
}
}
if(olbo$Coordinates == "EM"){
olb = olbo$estimObject
catOrdBiplotPenal = GetOrdinalBiplotObjectPenal(olbo$dataSet$ColumNames,olb,planex,planey)
}else{
if(olbo$Coordinates == "MIRT"){
rows = olbo$estimObject
catOrdBiplot = GetOrdinalBiplotObjectMIRT(rows,planex,planey)
}else{
stop("Coordinates for the items has not been specified.")
}
}
if(PlotVars){
if(olbo$Coordinates == "EM"){
D = 1
plot.ordinalBiplotPenal(catOrdBiplotPenal,olbo$NumFactors,D,planex,planey,xi=xlimi,xu=xlimu,
yi=ylimi,yu=ylimu,margin = margin, CexVar = CexVar,ColorVar = ColorVar,
PchVar = PchVar,levelsVar = olbo$dataSet$LevelNames)
}else{
D = 1.702
plot.ordinalBiplotPenal(catOrdBiplot,olbo$NumFactors,D,planex,planey,xi=xlimi,xu=xlimu,
yi=ylimi,yu=ylimu,margin = margin, CexVar = CexVar,ColorVar = ColorVar,
PchVar = PchVar,levelsVar = olbo$dataSet$LevelNames)
}
if(legendPlot){
legend("bottomright", legend=VarNames, col= ColorVar,pch=PchVar,cex=0.7)
}
}
if (PlotClus) {
if (is.null(Clusters)){
Clusters=as.factor(ones(c(n,1)))
}else{
Clusters=Clusters
}
A = cbind(olbo$RowCoords[,planex],olbo$RowCoords[,planey])
if(!is.null(colorCluster)){
col = colorCluster
}else{
col = colors()[20 + 2*c(1:length(levels(Clusters)))]
}
if(!is.null(ConfidentLevel)){
if((ConfidentLevel < 0) | (ConfidentLevel > 1)){
stop("The ConfidentLevel parameter has not a valid value. It should be between 0 and 1")
}
}
PlotClusters(A, Clusters, colors = col, chulls = chulls,centers = centers,ConfidentLevel=ConfidentLevel)
}
if(showIIC){
for(v in 1:ncol(olbo$dataSet$datanom)){
nameVariable = VarNames[v]
if(olbo$Coordinates == "EM"){
coeffic = catOrdBiplotPenal$matBiplot[,v]$coef
slopeort = catOrdBiplotPenal$matBiplot[,v]$slope
D = 1
}else{
if(olbo$Coordinates == "MIRT"){
coeffic = catOrdBiplot$matBiplot[,v]$coef
slopeort = catOrdBiplot$matBiplot[,v]$slope
D = 1.702
}else{
stop("Coordinates for the items has not been specified.")
}
}
numcat = max(olbo$dataSet$datanom[,v])
xi = iicxi
xu = iicxu
plotCurvesCategoriesVariable(coeffic,slopeort,D,numcat,nameVariable,xi,xu,planex,planey)
}
}
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.