2008年11月27日星期四

R中glm模型的实现(总结)

(1)因变量>20&<20
>data("plasma", package = "HSAUR")
>plasma_glm_1 <- glm(ESR ~ fibrinogen, data = plasma,family = binomial())
>plasma_glm_2 <- glm(ESR ~ fibrinogen + globulin, data = plasma,family = binomial())

画阴影图:
> data("plasma", package = "HSAUR")
> layout(matrix(1:2, ncol = 2))
> cdplot(ESR ~ fibrinogen, data = plasma)
> cdplot(ESR ~ globulin, data = plasma)
比较两个模型:anova(plasma_glm_1, plasma_glm_2, test = "Chisq")
比较1:confint(plasma_glm_1, parm = "fibrinogen")
2.exp(coef(plasma_glm_1)["fibrinogen"])
3.exp(confint(plasma_glm_1, parm = "fibrinogen"))

泡泡图:
> prob <- predict(plasma_glm_2, type = "response")
画点> plot(globulin ~ fibrinogen, data = plasma, xlim = c(2, 6),ylim = c(25, 55), pch = ".")
画圈> symbols(plasma$fibrinogen, plasma$globulin, circles = prob,add = TRUE)

(2)women's role ----fm1 <- cbind(agree, disagree) ~ sex + education
一共有两种情况,(1)是只做sex和education的主效应,(2)是主效应再加上两者的交互作用
> data("womensrole", package = "HSAUR")
> fm1 <- cbind(agree, disagree) ~ sex + education
> womensrole_glm_1 <- glm(fm1, data = womensrole,family = binomial())
> fm2 <- cbind(agree,disagree) ~ sex * education
> womensrole_glm_2 <- glm(fm2, data = womensrole,family = binomial())
*做预测曲线
(1)> role.fitted1 <- predict(womensrole_glm_1, type = "response")
(2)> role.fitted2 <- predict(womensrole_glm_2, type = "response")

先做一个语句>myplot <- function(role.fitted) {
f <- womensrole$sex == "Female"
plot(womensrole$education, role.fitted, type = "n",
ylab = "Probability of agreeing",
xlab = "Education", ylim = c(0,1))
lines(womensrole$education[!f], role.fitted[!f], lty = 1)
lines(womensrole$education[f], role.fitted[f], lty = 2)
lgtxt <- c("Fitted (Males)", "Fitted (Females)")
legend("topright", lgtxt, lty = 1:2, bty = "n")
y <- womensrole$agree / (womensrole$agree +
womensrole$disagree)
text(womensrole$education, y, ifelse(f, "\\VE", "\\MA"),
family = "HersheySerif", cex = 1.25)
}
再输出图形
(1)> myplot(role.fitted1)
(2)> myplot(role.fitted2)
*做残差散点图
> res <- residuals(womensrole_glm_2, type = "deviance")
> plot(predict(womensrole_glm_2), res,xlab="Fitted values", ylab = "Residuals",ylim = max(abs(res)) * c(-1,1))
> abline(h = 0, lty = 2)



(3)number ~ treat + age,
计数数据的两个问题:都是正值,不是正态分布。考虑用Poisson回归的GML,这种标准化线性模型,即对数函数来保证都是正值,且保证 Poisson error distribution另外的方法就是用quasi-likelihood
3.1Poisson回归
> data("polyps", package = "HSAUR")
> polyps_glm_1 <- glm(number ~ treat + age, data = polyps,family = poisson())
> summary(polyps_glm_1)
3.2quasi-likelihood回归
> data("polyps", package = "HSAUR")
> polyps_glm_2 <- glm(number ~ treat + age, data = polyps, family = quasipoisson())
> summary(polyps_glm_2)

**3.1比3.2的要更加显著,虽然两者都是显著的

2008年11月22日星期六

blogger中latex的使用

第一步:
在CSS中写
img.eq_latex {
padding: 0;
margin: 0;
border: 0;
}
在]]>之前,用来消除公式的边框
第二步:
参考http://doc.yourequations.com/,第一步是原文上没有的。
第三步:
可按latex为关键词找blogger中的插件,还不错


一个例子:
公式是
\inline \int_{-\infty}^{\infty}e^{-x^{2}}\;dx=\sqrt{\pi}
前面的代码是code lang="eq.latex",后面是/code,都要用<>括起来,再把公式放在中间
另外,\inline 语句可使插入美观
例子:
The equation using CODE tag with \inline command,
\inline \int_{-\infty}^{\infty}e^{-x^{2}}\;dx=\sqrt{\pi}
, and here is the following text.


\inline \sum_{\alpha}^{\beta}


第一个代码比较难找,找到了还不对,后来我把代码改了一下才能用。这样,效果也就不逊于WORDPRESS了
为虾米不能在google阅读器中显示那?who can help me??

2008年11月21日星期五

(转)李老师修改过的用SPSS做调节和中介的方法

转自沐辰http://10816038fmc.blogspot.com/~
先用descriptive statistics报告出M和SD,然后将所有的控制变量(人口学因素)、自变量和调节或中介变量以及因变量放入SPSS计算相关。得到相关表,发现与自变量和调节或中介或因变量相关显著的变量就要作为回归中需要考虑的变量放入回归方程中计算。
##这是预研究阶段的操作。在验证性研究中,不应当在正式研究的行文中出现这种data snoop的措辞。要时刻留心将模型的提出和模型的验证区别开来。如果未能实现这一点,意味着只得到提出模型的研究,没有得到验证模型的研究。我们在验证性研究中screen数据,目的是为了确保模型条件(异常值、分布)没有出现明显反常识的症状;而不是为了提出模型(设置哪些变量作为IV,是否设置交互项)。如果一套数据被用于提出模型,就不应当再被用于验证同一个模型。虽然国内很多学术刊物发表了类似的错误方法案例。

调节:自变量和调节变量要中心化:descriptive statistics –save standardized values as variables(勾选)
##这是标准化而不只是中心化。中心化是让变量减去自身的均值。标准化将让结果的解读丧失部分信息。
第一步放控制变量和调节变量(中心化以后的)
2 加入自变量(中心化后的)
3 加入自变量和调节变量的乘积(compute)(这些都是用回归中block那个选项卡,直接next,在里面放入多加入的变量就可以了)
##
##技术上的细节:中心化的目的是为了让截距项方便解读。仔细分析,课讨论是否必要对调节变量作中心化;如果自变量是nominal变量得到的dummy变量,是否需要中心化,

中介:不用中心化 ,同上面的1 2 步。

都是看R平方的改变和F的改变(sig)是否小于0.01或者0.05。
##中介效应是作sobel检验或者更精密的a*b置信区间,不是看R^2改变的显著性
##调节效应应当作相应的f^2置信区间,特别是在样本量大的时候。我们不能等待SPSS开发了这项功能再学习这项技术

画高低图的时候 把控制变量去掉重新做回归方程。用自变量和调节变量放入方程, 取(1,1)(-1,1)(1,-1)(-1,-1)这四个点画图。

##请讨论增改这份笔记。有三个方面可以增益的工作:
1. 加入操作图,或者syntax。二者有一样就有manual的价值
2. 操作上的错误或者模糊的地方订正
3. 联系挖掘更深一层的原理

目前的水准与国内同主题的网络资料没有差异,尚需进一步加工,才能达到预期的北大心理系研究生学习水准。
大家有任何增改意见都可以提出来哈~

(转)高级心理统计课考核方式通知

高级心理统计课考核方式通知:平时分值与笔试分值组成比例与基础心理统计课相同,但平时作品的形式更为灵活:
在内容上,可以写个案疑难分析、答疑学习笔记、在线问卷实战技巧总结等等,优劣主要由任课教师主观评判,主要考虑其他同学如果阅读该页面在学习上可能有多大受益。在这个变量被Covariate之后,篇幅不再有影响。
1.
用任何一种形式在网络上以学号公开发布,并将访问链接提交给学习委员。请学习委员将这些访问链接写入一份google
docs,该页面将公开发布。
2. 也可以用任何文件格式作为邮件附件提交给我。
其中,第一种方式包括但不限于:google
docs在线文档;个人blog。这两种形式的作品我会尽量通过在线协作评改或comment。第二种方式提交的作品,可以通过来办公室答疑获得评改。我也将在优选之后帮助发布在网页上。
与基础课不同的,高级课的作业中,任何引用部分都必须严格标出来源,严格区分创作与引用。适合提交给基础课的作业如果作为高级课作业提交,将不会减少也不会增益高级课的得分。重复提交给基础课作业与高级课作业的工作,将不会减少也不会增益高级课的得分。
我期待有1/10~1/5的同学有兴趣从这次课开始维护一个反映学术工作成果与学习心得笔记的长期个人博客。天道酬勤,时间+足迹才可以证明你们比所有稍微懒惰一点点的竞争者更为优秀。
同时,我也预期至少有1/3左右的同学无视本通知。
LI, Xiaoxu
我的VOA听写积分