PCA - scores

classic Classic list List threaded Threaded
10 messages Options
Reply | Threaded
Open this post in threaded view
|

PCA - scores

sclare
I am running a PCA, but would like to rotate my data and limit the  
number of factors that are analyzed.  I can do this using the  
"principal" command from the psych package [principal(my.data,  
nfactors=3,rotate="varimax")], but the issue is that this does not  
report scores for the Principal Components the way "princomp" does.

My question is:

Can you get an output of scores using "principal" OR, is there a way  
to limit the number of factors that are included when you use  
"princomp"?

Thanks,
Shari Clare

PhD Candidate
Department of Renewable Resources
University of Alberta
[hidden email]
780-492-2540








        [[alternative HTML version deleted]]

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: PCA - scores

Joshua Wiley-2
Hi Shari,

Yes, please look at the documentation for principal.  You can access
this (assuming you have loaded psych) by typing at the console:

?principal

note the logical argument "scores".

Here is a small example:

##############################
require(psych)
require(GPArotation)

dat <- principal(mtcars[, c("mpg", "hp", "wt")], nfactors = 1,
  rotate = "oblimin", scores = TRUE)

dat$scores
##############################

Cheerio,

Josh

On Thu, Mar 3, 2011 at 1:02 PM, Shari Clare <[hidden email]> wrote:

> I am running a PCA, but would like to rotate my data and limit the
> number of factors that are analyzed.  I can do this using the
> "principal" command from the psych package [principal(my.data,
> nfactors=3,rotate="varimax")], but the issue is that this does not
> report scores for the Principal Components the way "princomp" does.
>
> My question is:
>
> Can you get an output of scores using "principal" OR, is there a way
> to limit the number of factors that are included when you use
> "princomp"?
>
> Thanks,
> Shari Clare
>
> PhD Candidate
> Department of Renewable Resources
> University of Alberta
> [hidden email]
> 780-492-2540
>
>
>
>
>
>
>
>
>        [[alternative HTML version deleted]]
>
> ______________________________________________
> [hidden email] mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



--
Joshua Wiley
Ph.D. Student, Health Psychology
University of California, Los Angeles
http://www.joshuawiley.com/

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: PCA - scores

William Revelle
Shari,
   Josh partly answered your question, but his example did not include
rotation because he took out just one factor.

Try:

require(psych)
mt.pc <- principal(mtcars,3,scores=TRUE)   #this gives you the
varimax rotated first 3 principal components
#pc.scores <- mt.pc$scores     #here are the scores

biplot(mt.pc)    #show the data as well as the principal components in a biplot



Bill


At 5:15 PM -0800 3/3/11, Joshua Wiley wrote:

>Hi Shari,
>
>Yes, please look at the documentation for principal.  You can access
>this (assuming you have loaded psych) by typing at the console:
>
>?principal
>
>note the logical argument "scores".
>
>Here is a small example:
>
>##############################
>require(psych)
>require(GPArotation)
>
>dat <- principal(mtcars[, c("mpg", "hp", "wt")], nfactors = 1,
>   rotate = "oblimin", scores = TRUE)
>
>dat$scores
>##############################
>
>Cheerio,
>
>Josh
>
>On Thu, Mar 3, 2011 at 1:02 PM, Shari Clare <[hidden email]> wrote:
>>  I am running a PCA, but would like to rotate my data and limit the
>>  number of factors that are analyzed.  I can do this using the
>>  "principal" command from the psych package [principal(my.data,
>>  nfactors=3,rotate="varimax")], but the issue is that this does not
>>  report scores for the Principal Components the way "princomp" does.
>>
>>  My question is:
>>
>>  Can you get an output of scores using "principal" OR, is there a way
>>  to limit the number of factors that are included when you use
>>  "princomp"?
>>
>>  Thanks,
>>  Shari Clare
>>
>>  PhD Candidate
>>  Department of Renewable Resources
>>  University of Alberta
>>  [hidden email]
>>  780-492-2540
>>
>>
>>
>>
>>
>>
>>
>>
>>         [[alternative HTML version deleted]]
>>
>>  ______________________________________________
>>  [hidden email] mailing list
>>  https://stat.ethz.ch/mailman/listinfo/r-help
>>  PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>>  and provide commented, minimal, self-contained, reproducible code.
>>
>
>
>
>--
>Joshua Wiley
>Ph.D. Student, Health Psychology
>University of California, Los Angeles
>http://www.joshuawiley.com/
>
>______________________________________________
>[hidden email] mailing list
>https://stat.ethz.ch/mailman/listinfo/r-help
>PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>and provide commented, minimal, self-contained, reproducible code.

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: PCA - scores

sclare
Hi Bill and Josh:

When I run any "principal" code with scores=TRUE, I get the following  
Error:

Error in principal (my.data,3,scores=TRUE) : unused argument  
(scores=TRUE)

Thoughts?

Thanks,
Shari






On 3-Mar-11, at 9:42 PM, William Revelle wrote:

> Shari,
>  Josh partly answered your question, but his example did not include  
> rotation because he took out just one factor.
>
> Try:
>
> require(psych)
> mt.pc <- principal(mtcars,3,scores=TRUE)   #this gives you the  
> varimax rotated first 3 principal components
> #pc.scores <- mt.pc$scores     #here are the scores
>
> biplot(mt.pc)    #show the data as well as the principal components  
> in a biplot
>
>
>
> Bill
>
>
> At 5:15 PM -0800 3/3/11, Joshua Wiley wrote:
>> Hi Shari,
>>
>> Yes, please look at the documentation for principal.  You can access
>> this (assuming you have loaded psych) by typing at the console:
>>
>> ?principal
>>
>> note the logical argument "scores".
>>
>> Here is a small example:
>>
>> ##############################
>> require(psych)
>> require(GPArotation)
>>
>> dat <- principal(mtcars[, c("mpg", "hp", "wt")], nfactors = 1,
>>  rotate = "oblimin", scores = TRUE)
>>
>> dat$scores
>> ##############################
>>
>> Cheerio,
>>
>> Josh
>>
>> On Thu, Mar 3, 2011 at 1:02 PM, Shari Clare <[hidden email]>  
>> wrote:
>>> I am running a PCA, but would like to rotate my data and limit the
>>> number of factors that are analyzed.  I can do this using the
>>> "principal" command from the psych package [principal(my.data,
>>> nfactors=3,rotate="varimax")], but the issue is that this does not
>>> report scores for the Principal Components the way "princomp" does.
>>>
>>> My question is:
>>>
>>> Can you get an output of scores using "principal" OR, is there a way
>>> to limit the number of factors that are included when you use
>>> "princomp"?
>>>
>>> Thanks,
>>> Shari Clare
>>>
>>> PhD Candidate
>>> Department of Renewable Resources
>>> University of Alberta
>>> [hidden email]
>>> 780-492-2540
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>        [[alternative HTML version deleted]]
>>>
>>> ______________________________________________
>>> [hidden email] mailing list
>>> https://stat.ethz.ch/mailman/listinfo/r-help
>>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>>> and provide commented, minimal, self-contained, reproducible code.
>>>
>>
>>
>>
>> --
>> Joshua Wiley
>> Ph.D. Student, Health Psychology
>> University of California, Los Angeles
>> http://www.joshuawiley.com/
>>
>> ______________________________________________
>> [hidden email] mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>


        [[alternative HTML version deleted]]

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: PCA - scores

William Revelle
At 9:52 AM -0700 3/4/11, Shari Clare wrote:
>Hi Bill and Josh:
>
>When I run any "principal" code with scores=TRUE, I get the following Error:
>
>Error in principal (my.data,3,scores=TRUE) : unused argument (scores=TRUE)
>
>
>Thoughts?

What version of psych are you using?

Does it work on the example I sent (see below)?



>
>Thanks,
>Shari
>
>
>
>
>
>
>On 3-Mar-11, at 9:42 PM, William Revelle wrote:
>
>>Shari,
>>  Josh partly answered your question, but his example did not
>>include rotation because he took out just one factor.
>>
>>Try:
>>
>>require(psych)
>>mt.pc <- principal(mtcars,3,scores=TRUE)   #this gives you the
>>varimax rotated first 3 principal components
>>#pc.scores <- mt.pc$scores     #here are the scores
>>
>>biplot(mt.pc)    #show the data as well as the principal components
>>in a biplot
>>
>>
>>
>>Bill
>>
>>
>>At 5:15 PM -0800 3/3/11, Joshua Wiley wrote:
>>
>>>Hi Shari,
>>>
>>>
>>>Yes, please look at the documentation for principal.  You can access
>>>
>>>this (assuming you have loaded psych) by typing at the console:
>>>
>>>
>>>?principal
>>>
>>>
>>>note the logical argument "scores".
>>>
>>>
>>>Here is a small example:
>>>
>>>
>>>##############################
>>>
>>>require(psych)
>>>
>>>require(GPArotation)
>>>
>>>
>>>dat <- principal(mtcars[, c("mpg", "hp", "wt")], nfactors = 1,
>>>
>>>  rotate = "oblimin", scores = TRUE)
>>>
>>>
>>>dat$scores
>>>
>>>##############################
>>>
>>>
>>>Cheerio,
>>>
>>>
>>>Josh
>>>
>>>
>>>On Thu, Mar 3, 2011 at 1:02 PM, Shari Clare
>>><<mailto:[hidden email]>[hidden email]> wrote:
>>>
>>>>I am running a PCA, but would like to rotate my data and limit the
>>>>
>>>>number of factors that are analyzed.  I can do this using the
>>>>
>>>>"principal" command from the psych package [principal(my.data,
>>>>
>>>>nfactors=3,rotate="varimax")], but the issue is that this does not
>>>>
>>>>report scores for the Principal Components the way "princomp" does.
>>>>
>>>>
>>>>My question is:
>>>>
>>>>
>>>>Can you get an output of scores using "principal" OR, is there a way
>>>>
>>>>to limit the number of factors that are included when you use
>>>>
>>>>"princomp"?
>>>>
>>>>
>>>>Thanks,
>>>>
>>>>Shari Clare
>>>>
>>>>
>>>>PhD Candidate
>>>>
>>>>Department of Renewable Resources
>>>>
>>>>University of Alberta
>>>>
>>>><mailto:[hidden email]>[hidden email]
>>>>
>>>>780-492-2540
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>        [[alternative HTML version deleted]]
>>>>
>>>>
>>>>______________________________________________
>>>>
>>>><mailto:[hidden email]>[hidden email] mailing list
>>>>
>>>><https://stat.ethz.ch/mailman/listinfo/r-help>https://stat.ethz.ch/mailman/listinfo/r-help
>>>>
>>>>PLEASE do read the posting guide
>>>><http://www.R-project.org/posting-guide.html>http://www.R-project.org/posting-guide.html
>>>>
>>>>and provide commented, minimal, self-contained, reproducible code.
>>>>
>>>>
>>>
>>>
>>>
>>>--
>>>
>>>Joshua Wiley
>>>
>>>Ph.D. Student, Health Psychology
>>>
>>>University of California, Los Angeles
>>>
>>><http://www.joshuawiley.com/>http://www.joshuawiley.com/
>>>
>>>
>>>______________________________________________
>>>
>>><mailto:[hidden email]>[hidden email] mailing list
>>>
>>><https://stat.ethz.ch/mailman/listinfo/r-help>https://stat.ethz.ch/mailman/listinfo/r-help
>>>
>>>PLEASE do read the posting guide
>>><http://www.R-project.org/posting-guide.html>http://www.R-project.org/posting-guide.html
>>>
>>>and provide commented, minimal, self-contained, reproducible code.

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: PCA - scores

Uwe Ligges-3
In reply to this post by sclare


On 04.03.2011 17:52, Shari Clare wrote:
> Hi Bill and Josh:
>
> When I run any "principal" code with scores=TRUE, I get the following
> Error:
>
> Error in principal (my.data,3,scores=TRUE) : unused argument
> (scores=TRUE)
>
> Thoughts?

Your psych version (and probably also your R version) is outdated?
Please upgrade both R and your packages.

Best,
Uwe Ligges




>
> Thanks,
> Shari
>
>
>
>
>
>
> On 3-Mar-11, at 9:42 PM, William Revelle wrote:
>
>> Shari,
>>   Josh partly answered your question, but his example did not include
>> rotation because he took out just one factor.
>>
>> Try:
>>
>> require(psych)
>> mt.pc<- principal(mtcars,3,scores=TRUE)   #this gives you the
>> varimax rotated first 3 principal components
>> #pc.scores<- mt.pc$scores     #here are the scores
>>
>> biplot(mt.pc)    #show the data as well as the principal components
>> in a biplot
>>
>>
>>
>> Bill
>>
>>
>> At 5:15 PM -0800 3/3/11, Joshua Wiley wrote:
>>> Hi Shari,
>>>
>>> Yes, please look at the documentation for principal.  You can access
>>> this (assuming you have loaded psych) by typing at the console:
>>>
>>> ?principal
>>>
>>> note the logical argument "scores".
>>>
>>> Here is a small example:
>>>
>>> ##############################
>>> require(psych)
>>> require(GPArotation)
>>>
>>> dat<- principal(mtcars[, c("mpg", "hp", "wt")], nfactors = 1,
>>>   rotate = "oblimin", scores = TRUE)
>>>
>>> dat$scores
>>> ##############################
>>>
>>> Cheerio,
>>>
>>> Josh
>>>
>>> On Thu, Mar 3, 2011 at 1:02 PM, Shari Clare<[hidden email]>
>>> wrote:
>>>> I am running a PCA, but would like to rotate my data and limit the
>>>> number of factors that are analyzed.  I can do this using the
>>>> "principal" command from the psych package [principal(my.data,
>>>> nfactors=3,rotate="varimax")], but the issue is that this does not
>>>> report scores for the Principal Components the way "princomp" does.
>>>>
>>>> My question is:
>>>>
>>>> Can you get an output of scores using "principal" OR, is there a way
>>>> to limit the number of factors that are included when you use
>>>> "princomp"?
>>>>
>>>> Thanks,
>>>> Shari Clare
>>>>
>>>> PhD Candidate
>>>> Department of Renewable Resources
>>>> University of Alberta
>>>> [hidden email]
>>>> 780-492-2540
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>         [[alternative HTML version deleted]]
>>>>
>>>> ______________________________________________
>>>> [hidden email] mailing list
>>>> https://stat.ethz.ch/mailman/listinfo/r-help
>>>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>>>> and provide commented, minimal, self-contained, reproducible code.
>>>>
>>>
>>>
>>>
>>> --
>>> Joshua Wiley
>>> Ph.D. Student, Health Psychology
>>> University of California, Los Angeles
>>> http://www.joshuawiley.com/
>>>
>>> ______________________________________________
>>> [hidden email] mailing list
>>> https://stat.ethz.ch/mailman/listinfo/r-help
>>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>>> and provide commented, minimal, self-contained, reproducible code.
>>
>
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> [hidden email] mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Color-coded Biplot of Varimax rotated factors from PCA based factor analysis

Barbara Doll
In reply to this post by William Revelle
I'm using R for PCA and? factor analysis. I want to create biplots of
varimax rotated factors that color-code points by their classification.
My research is on streams that are urban and rural. So, I want to color
code them by this classification. If you just do a biplot from prcomp or
princomp, you cannot add this color. So, I have used some code developed
by a graduate student in our statistics department here at NC State
University. However, when you compare the two biplots, the observed
points are not in the same location. The variable vectors match up, but
not the points. I'm not sure why. The code is below. Please help.


label=data[,"Urban.Rural"]
indexU<-which(label=="U")
indexR<-which(label=="R")

collab<-rep(0,length(data[,1]))
collab[indexU]<-"Blue"
collab[indexR]<-"Green"

library(psych)

fit <- principal(mydata, nfactors = num.fac, rotate="varimax", scores =
TRUE)

z1 <- sum(fit2$loadings[,1]^2)  ### need to scale scores and loadings by
these factors
z2 <- sum(fit2$loadings[,2]^2)

biplot(fit$scores[,c(1,2)]/c(z1, z2), loadings(fit)[,c(1,2)]*c(z1, z2),
xlabs=rep("", length(collab)),  col = c("black", "orange"))
legend(x="bottomright", legend=c("Urban","Rural"),
text.col=c("Blue","Green"), bg="Grey90")

### scale the plot parameters
rrr<-apply(fit$scores[,1:2],2, range)
(abs(rrr)+.1)*sign(rrr)
par(usr=as.vector(rrr))

### now include the colored points
points(fit$scores[,c(1,2)], col=collab, pch=20)




        [[alternative HTML version deleted]]

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: Color-coded Biplot of Varimax rotated factors from PCA based factor analysis

PIKAL Petr
Hi

Without some data it is difficult to say where is difference. However from code below it seems to me that it is an adaptation from my attempt to produce colour coded points on biplot.

http://tolstoy.newcastle.edu.au/R/e4/help/08/08/20012.html

fit<-princomp(iris[,1:4], cor=T)
biplot(fit, xlabs=rep("", nrow(iris)))
rrr<-apply(fit$scores[,1:2],2, range)
par(usr=as.vector(rrr))
points(fit$scores[,1:2], col=rainbow(4)[(as.numeric(iris[,5]))], pch=20)

This results in presenting points coloured according to some rules but the points are slightly shifted from those made by biplot itself. So AFAIK you either can use biplot and no colouring or use the code to present faked coloured points slightly shifted.

Compare with just

biplot(fit)
points(fit$scores[,1:2], col=rainbow(4)[(as.numeric(iris[,5]))], pch=20)

But maybe something changed from 2008 and now biplot with coloured point is available somewhere.

Regards
Petr



>
> I'm using R for PCA and? factor analysis. I want to create biplots of
> varimax rotated factors that color-code points by their classification.
> My research is on streams that are urban and rural. So, I want to color
> code them by this classification. If you just do a biplot from prcomp
> or princomp, you cannot add this color. So, I have used some code
> developed by a graduate student in our statistics department here at NC
> State University. However, when you compare the two biplots, the
> observed points are not in the same location. The variable vectors
> match up, but not the points. I'm not sure why. The code is below.
> Please help.
>
>
> label=data[,"Urban.Rural"]
> indexU<-which(label=="U")
> indexR<-which(label=="R")
>
> collab<-rep(0,length(data[,1]))
> collab[indexU]<-"Blue"
> collab[indexR]<-"Green"
>
> library(psych)
>
> fit <- principal(mydata, nfactors = num.fac, rotate="varimax", scores =
> TRUE)
>
> z1 <- sum(fit2$loadings[,1]^2)  ### need to scale scores and loadings
> by these factors
> z2 <- sum(fit2$loadings[,2]^2)
>
> biplot(fit$scores[,c(1,2)]/c(z1, z2), loadings(fit)[,c(1,2)]*c(z1, z2),
> xlabs=rep("", length(collab)),  col = c("black", "orange"))
> legend(x="bottomright", legend=c("Urban","Rural"),
> text.col=c("Blue","Green"), bg="Grey90")
>
> ### scale the plot parameters
> rrr<-apply(fit$scores[,1:2],2, range)
> (abs(rrr)+.1)*sign(rrr)
> par(usr=as.vector(rrr))
>
> ### now include the colored points
> points(fit$scores[,c(1,2)], col=collab, pch=20)
>
>
>
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> [hidden email] mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-
> guide.html
> and provide commented, minimal, self-contained, reproducible code.

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: Color-coded Biplot of Varimax rotated factors from PCA based factor analysis

Adams, Jean
In reply to this post by Barbara Doll
I suggest that you re-post your question, but this time make it easy for
readers to help you by providing a simple reproducible example of your
problem.  Provide some fake data, or a small subset of your data using,
for example, the function dput().  Trim the fat from your code.  Get rid
of any arguments or functions that are not necessary to demonstrate your
problem.  Include in your code the function library() or require() to load
any needed packages.  If you make it so that a reader can simply copy and
paste your code into R, you will get more help.

Jean


Barbara Doll <[hidden email]> wrote on 08/15/2012 01:29:09 PM:
>
> I'm using R for PCA and? factor analysis. I want to create biplots of
> varimax rotated factors that color-code points by their classification.
> My research is on streams that are urban and rural. So, I want to color
> code them by this classification. If you just do a biplot from prcomp or

> princomp, you cannot add this color. So, I have used some code developed

> by a graduate student in our statistics department here at NC State
> University. However, when you compare the two biplots, the observed
> points are not in the same location. The variable vectors match up, but
> not the points. I'm not sure why. The code is below. Please help.
>
>
> label=data[,"Urban.Rural"]
> indexU<-which(label=="U")
> indexR<-which(label=="R")
>
> collab<-rep(0,length(data[,1]))
> collab[indexU]<-"Blue"
> collab[indexR]<-"Green"
>
> library(psych)
>
> fit <- principal(mydata, nfactors = num.fac, rotate="varimax", scores =
> TRUE)
>
> z1 <- sum(fit2$loadings[,1]^2)  ### need to scale scores and loadings by

> these factors
> z2 <- sum(fit2$loadings[,2]^2)
>
> biplot(fit$scores[,c(1,2)]/c(z1, z2), loadings(fit)[,c(1,2)]*c(z1, z2),
> xlabs=rep("", length(collab)),  col = c("black", "orange"))
> legend(x="bottomright", legend=c("Urban","Rural"),
> text.col=c("Blue","Green"), bg="Grey90")
>
> ### scale the plot parameters
> rrr<-apply(fit$scores[,1:2],2, range)
> (abs(rrr)+.1)*sign(rrr)
> par(usr=as.vector(rrr))
>
> ### now include the colored points
> points(fit$scores[,c(1,2)], col=collab, pch=20)

        [[alternative HTML version deleted]]

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

using prcomp() and varimax()

MikeAmato.WI
This post has NOT been accepted by the mailing list yet.
Hello,
I am attempting to do a principal components analysis on 15 survey items. I want to use a varimax rotation on the retained components, but I am dubious of the output I am getting, and so I suspect I am doing something wrong. I proceed in the following steps:

 1) use prcomp() to inspect all 15 components, and decide which to retain
 2) run prcomp() again, using the "tol" parameter to omit unwanted components
 3) pass the output of step 2 to varimax()

My concern is with the reported proportions of variance for the 3 components after varimax rotation. It looks like each of my 3 components explains 1/15 of the total variance, summing to a cumulative proportion of 20% of variance explained. But those 3 components I retained should now be the only components in the analysis, so they should be able to account for 100% of the explained variance.

I am able to get reliable seeming results using principal() from the "psych" package, in which the total amount of variance explained by my retained components does not differ before or after rotation. But principal() uses varimax(), so I suspect I am either doing something wrong or misinterpreting the output when using the base package functions.  

Am I doing something wrong when attempting to retain only 3 components?
Am I using varimax() incorrectly?
Am I misinterpreting the returned values from varimax()?

Thanks for any help,
Mike



Here is a link to the data file I am using: https://www.dropbox.com/s/scypebzy0nnhlwk/pca_sampledata.txt

### step 1 ###
> d1 = read.table("pca_sampledata.txt", T)
> m1 = with(d1, ~ g.enjoy + g.look + g.cost + g.fit + g.health + g.resale + b.withstand + b.satisfy + b.vegetated + b.everyone + b.harmed + b.eco + b.ingenuity + b.security + b.proud)
> pca1 = prcomp(m1)
> summary(pca1) #output truncated for this posting
Importance of components:
                          PC1    PC2    PC3     PC4     PC5 ...    PC15
Standard deviation     1.5531 1.3064 1.1695 0.93512 0.92167 ... 0.35500
Proportion of Variance 0.2199 0.1556 0.1247 0.07972 0.07744 ... 0.01149
Cumulative Proportion  0.2199 0.3755 0.5002 0.57988 0.65732 ... 1.00000


### step 2 ###
> pca2 = prcomp(m1, tol=.75)
> summary(pca2) #full output shown
Importance of components:
                          PC1    PC2    PC3
Standard deviation     1.5531 1.3064 1.1695
Proportion of Variance 0.4397 0.3111 0.2493
Cumulative Proportion  0.4397 0.7507 1.0000


### step 3 ###
> pca3 = varimax(pca2$rotation)
> pca3
> ...
>                  PC1   PC2   PC3
> SS loadings    1.000 1.000 1.000
> Proportion Var 0.067 0.067 0.067
> Cumulative Var 0.067 0.133 0.200