现在的位置: 首页Skills 点滴>正文
数量化理论1程序-R语言
2011年04月03日 Skills 点滴 暂无评论

R语言下载地址

 

corpcor程序包 自行下载加载(通过界面上的 程序包)

复制下面文本到R语言控制台(console)中 粘贴 回车 可以看到示例结果

以下内容转自普蘭塔论坛 作者gagea

  1. #数量化方法I
  2. library(corpcor)
  3. y=c(3,5,6,7,9,11,9,7,7,6)
  4. c11=c(1,0,0,1,0,0,1,0,0,0)
  5. c12=c(0,1,0,0,1,0,0,1,0,1)
  6. c13=c(0,0,1,0,0,1,0,0,1,0)
  7. c21=c(1,1,1,0,0,0,0,0,1,1)
  8. c22=c(0,0,0,1,1,1,1,1,0,0)
  9. test=data.frame(y=y,c11=c11,c12=c12,c13=c13,c21=c21,c22=c22)
  10. #定义数量化理论I分析函数
  11. ##x为反应矩阵,用数据框表示;y为基准变量,用向量表示;nItem为项目数;nCati为各项目类目数,用向量表示。
  12. ##ItemNames为项目名称向量。
  13. TOQO=function(x,y,nItem,ItemNames,nCati){
  14. xnames=names(x)
  15. xm=as.matrix(x)
  16. ym=as.matrix(y)
  17. x.x=t(xm)%*%xm#正则方程系数矩阵
  18. x.y=t(xm)%*%ym#正则方程右端常数向量
  19. #对每个j=2,...,m,删去第j项目第一类目的方程
  20. Index=1:sum(nCati)
  21. IndexDel=NULL
  22. for(i in 2:nItem) {Index=Index[Index!=sum(nCati[1:(i-1)])+1]#确定保留索引号
  23. IndexDel[i-1]=sum(nCati[1:(i-1)])+1#确定要删除的索引号
  24. }
  25. x.x_2=x.x[Index,Index]
  26. x.y_2=x.y[Index]
  27. b1=solve(x.x_2,x.y_2)#预测方程部分解
  28. b=NULL#初始化预测方程完整解。
  29. b[Index]=b1
  30. b[IndexDel]=0
  31. ItemM=matrix(nrow=nrow(x),ncol=nItem)#项目预测得分矩阵
  32. for(i in 1:nItem){
  33. if(i==1) ItemM[,i]=xm[,1:nCati[i]]%*%b[1:nCati[i]]
  34. else ItemM[,i]=xm[,IndexDel[i-1]:(IndexDel[i-1]+nCati[i]-1)]%*%b[IndexDel[i-1]:(IndexDel[i-1]+nCati[i]-1)]
  35. }
  36. M=cbind(ym,ItemM)
  37. attributes(M)$dimnames[[2]]=c("y",ItemNames)
  38. xcor=cor(M)#相关系数矩阵
  39. xpcor=cor2pcor(xcor)#偏相关矩阵
  40. pcor=xpcor[2:(nItem+1),1]#偏相关系数
  41. pcordf=as.data.frame(t(pcor))
  42. ItemVar=apply(ItemM,2,var)
  43. yVar=var(ym)
  44. ItemVarRatio=ItemVar/yVar#方差比
  45. IVR=as.data.frame(t(ItemVarRatio))
  46. Range=NULL
  47. for(i in 1:nItem){#范围
  48. if(i==1) Range[i]=max(b[1:nCati[i]])-min(b[1:nCati[i]])
  49. else Range[i]=max(b[IndexDel[i-1]:(IndexDel[i-1]+nCati[i]-1)])-min(b[IndexDel[i-1]:(IndexDel[i-1]+nCati[i]-1)])
  50. }
  51. Range=as.data.frame(t(Range))
  52. Table=rbind(pcordf,IVR,Range)
  53. names(Table)=ItemNames
  54. row.names(Table)=c("偏相关系数","方差比","范围")
  55. ans=list(b=b,Table=Table)
  56. ans
  57. }
  58. #例:
  59. x=test[,2:6]
  60. y=test[,1]
  61. nItem=2
  62. nCati=c(3,2)
  63. ItemNames=c("体重","性别")
  64. TOQO(x,y,nItem,ItemNames,nCati)
0

给我留言

留言无头像?

×
腾讯微博