Dear all,
I'm writing my manuscript to publish after analysis my final data with
ANOVA, ANCOVA, MANCOVA. In a section of my result, I did correlation of my
data (2 categirical factors with 2 levels: Quantity & Quality; 2 dependent
var: Irid.area & Casa.PC1, and 1 co-var: SL). But as some traits (here
Irid.area) are significantly influenced by the covariate (standard length,
SL), I need to use the partial correlation. I know how to calculate it with
JMP, but as I used R to analyse all of my data (first time in my life) for
this manuscript, can anyone help me to find a solution for this problem? I
got some libraries to calculate it (e.g. ppcor, ggm, etc.), but none of
them fit to my required analysis (fitting covariate and subset group) in
the model.
Any help will be very much appreciated.
> ### Datafrmame> data1 <- read.csv(file.choose(),header=TRUE)
#Partial correlation.csv> data1 Quantity Quality SL Irid.area
Casa.PC1
1 High Low 16.38 10.31 1.711739555
2 High High 15.95 16.52 0.013383537
3 High High 15.69 12.74 2.228490878
4 High Low 14.76 9.80 1.554975833
5 High Low 14.63 12.95 1.823767970
6 High High 14.32 14.21 3.152059841
7 High High 14.95 12.57 2.069265040
8 High Low 15.37 13.55 1.886027422
9 High Low 14.73 14.18 1.127440602
10 High High 16.08 15.98 1.435563307
11 High High 15.78 16.76 2.433261686
12 High Low 15.22 12.12 0.927454986
13 High Low 14.22 10.91 2.328899576
14 High High 14.47 11.03 1.522923487
15 High Low 13.98 10.03 2.342535074
16 High Low 14.99 11.44 0.749529924
17 High High 16.51 20.16 2.993905677
18 High High 14.83 16.82 2.227315597
19 High Low 15.17 19.21 1.685063793
20 High Low 16.29 20.31 1.551704440
21 High High 16.23 15.03 1.982319336
22 High High 14.18 14.80 1.839910851
23 High Low 16.11 12.92 1.443240647
24 High Low 13.95 7.60 2.034192171
25 High High 17.54 17.80 2.188306237
26 High Low 16.24 19.29 1.531264746
27 High High 14.79 12.98 1.465644134
28 High Low 15.87 14.85 1.372494892
29 High High 16.09 13.71 1.462037152
30 High Low 14.34 13.53 1.365588960
31 High High 14.93 12.91 0.729212386
32 High High 15.89 16.98 0.136175317
33 High Low 16.11 11.93 1.442761666
34 High Low 15.25 15.49 0.834442777
35 High High 15.84 17.65 1.471713978
36 High High 15.61 18.00 1.949457500
37 High Low 15.42 13.87 0.200098471
38 High Low 14.91 11.23 0.981988071
39 High High 15.69 5.74 -0.445941360
40 High High 15.13 9.07 1.387947896
41 High Low 15.04 15.87 1.480980400
42 High Low 17.08 17.24 2.620029423
43 High High 15.85 12.47 0.027278890
44 High High 15.35 10.44 2.597373230
45 High Low 15.62 12.11 0.030653396
46 High High 17.96 17.50 1.544922124
47 High Low 17.25 17.87 1.705053951
48 High Low 15.56 19.72 1.688867665
49 High High 16.27 13.15 0.111371757
50 High Low 16.68 15.43 1.538012366
51 High High 15.78 15.07 0.744555741
52 Low High 14.72 13.34 -0.682505420
53 Low Low 14.93 14.07 -1.641494605
54 Low High 13.94 13.22 -1.172268647
55 Low High 14.01 18.65 -0.996656064
56 Low Low 14.33 17.16 -1.789728167
57 Low Low 14.57 12.43 -0.827526343
58 Low High 14.01 15.29 -1.350691602
59 Low Low 14.22 16.98 -1.688278221
60 Low High 13.45 14.40 -1.182117327
61 Low High 13.44 16.57 -1.358976542
62 Low Low 14.76 15.58 0.334534454
63 Low Low 14.85 17.65 0.251766383
64 Low High 13.42 10.99 -0.526634460
65 Low High 14.07 16.88 -1.112579922
66 Low Low 14.15 16.41 -0.971918177
67 Low Low 14.78 11.95 -1.179074800
68 Low High 14.84 17.62 -0.777057705
69 Low High 15.16 14.09 -1.224388816
70 Low Low 14.60 15.03 -0.775478528
71 Low High 13.74 10.01 -0.917153842
72 Low High 13.54 12.34 -0.822895877
73 Low Low 14.04 11.86 0.002789116
74 Low High 15.73 18.50 -1.209469875
75 Low Low 15.14 16.85 -0.479090055
76 Low Low 14.86 17.32 -1.897204235
77 Low High 14.43 11.20 0.469569392
78 Low Low 14.01 15.55 -1.025059269
79 Low High 14.20 11.67 -0.770451072
80 Low High 16.16 17.34 -0.274527631
81 Low Low 14.63 13.52 -1.070187945
82 Low Low 15.83 14.85 -1.627211162
83 Low High 14.70 14.81 -1.694118608
84 Low High 13.91 14.48 -1.635459183
85 Low Low 13.95 16.05 -1.449612666
86 Low Low 14.03 12.58 -1.685968841
87 Low High 14.82 13.57 -0.097426417
88 Low High 14.32 12.16 -1.403512009
89 Low Low 14.33 7.66 -1.336654713
90 Low Low 15.01 10.15 -1.257019268
91 Low High 14.01 9.79 -0.715404495
92 Low Low 14.25 17.38 -1.296954022
93 Low High 14.55 16.11 -0.616895943
94 Low High 13.98 11.49 -0.654017365
95 Low Low 15.59 8.43 -1.708330027
96 Low Low 15.02 16.88 -1.352913634
97 Low High 13.99 9.64 -0.499793618
98 Low Low 13.98 12.25 -1.265336955
99 Low High 13.94 13.79 -0.263925513
100 Low Low 15.03 20.39 -0.720121308
101 Low Low 13.93 14.63 -0.908570400> ### COrrelation
test according to group> library("MASS")> with(data1, cor.test(~
Irid.area + Casa.PC1, subset=(Quantity=="High")))# gives cor, df+2,
p-values
Pearson's product-moment correlation
data: Irid.area and Casa.PC1
t = 1.5795, df = 49, p-value = 0.1206
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
-0.05905155 0.46734855
sample estimates:
cor
0.2201142> with(data1, cor.test(~ Irid.area + Casa.PC1,
subset=(Quantity=="Low")))# gives cor, df+2, p-values
Pearson's product-moment correlation
data: Irid.area and Casa.PC1
t = -0.4275, df = 48, p-value = 0.6709
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
-0.3342116 0.2205349
sample estimates:
cor
-0.06159377> #### Effect size from two-way ANOVA ####> anova<- aov(Irid.area ~
Quantity*Quality+SL, data=data1)> summary(anova) Df Sum Sq
Mean Sq F value Pr(>F)
Quantity 1 0.0 0.04 0.004 0.947
Quality 1 0.3 0.26 0.032 0.859
SL 1 149.5 149.49 18.027 5.03e-05 ***
Quantity:Quality 1 0.2 0.18 0.022 0.883
Residuals 96 796.1 8.29
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 >
etaSquared( anova ) # effect size eta.sq
eta.sq.part
Quantity 5.682119e-02 0.0632542041
Quality 6.577118e-05 0.0000781554
SL 1.510031e-01 0.1521470868
Quantity:Quality 1.922552e-04 0.0002284211> ### partial correlation
(pcor) tests:> library(ggm)Loading required package: graphError in
loadNamespace(i[[1L]], c(lib.loc, .libPaths())) :
there is no package called ‘BiocGenerics’In addition: Warning
messages:1: package ‘ggm’ was built under R version 2.15.3 2: package
‘graph’ was built under R version 3.0.1 Error: package ‘graph’ could
not be loaded> data2<- data1[, c("Irid.area",
"Casa.PC1", "SL"),
Quantity == "High"]> pcor(c("Irid.area",
"Casa.PC1",
"SL"),var(data2))Error in match.arg(method) : 'arg' must be
NULL or a
character vector> pc<-pcor(data2)> pc$estimate
Irid.area Casa.PC1 SL
Irid.area 1.0000000 -0.1313475 0.3387663
Casa.PC1 -0.1313475 1.0000000 0.5061438
SL 0.3387663 0.5061438 1.0000000
$p.value
Irid.area Casa.PC1 SL
Irid.area 0.0000000000 1.896426e-01 3.647247e-04
Casa.PC1 0.1896425857 0.000000e+00 6.258573e-09
SL 0.0003647247 6.258573e-09 0.000000e+00
$statistic
Irid.area Casa.PC1 SL
Irid.area 0.000000 -1.311637 3.564375
Casa.PC1 -1.311637 0.000000 5.809698
SL 3.564375 5.809698 0.000000
$n
[1] 101
$gp
[1] 1
$method
[1] "pearson"
Cheers,
--
MD. MOSHIUR RAHMAN
PhD Candidate
School of Animal Biology/Zoology (M092)
University of Western Australia
35 Stirling Hwy, Crawley, WA, 6009
Australia.
Mob.: 061-425205507
[[alternative HTML version deleted]]
Please don't cross post [1].
Please post minimal example data only. You should not be trying to (or appearing
to) use us to get your work done.
[1] http://www.r-project.org/mail.html
---------------------------------------------------------------------------
Jeff Newmiller The ..... ..... Go Live...
DCN:<jdnewmil at dcn.davis.ca.us> Basics: ##.#. ##.#. Live
Go...
Live: OO#.. Dead: OO#.. Playing
Research Engineer (Solar/Batteries O.O#. #.O#. with
/Software/Embedded Controllers) .OO#. .OO#. rocks...1k
---------------------------------------------------------------------------
Sent from my phone. Please excuse my brevity.
Moshiur Rahman <mrahmankufmrt at gmail.com> wrote:>Dear all,
>
>I'm writing my manuscript to publish after analysis my final data with
>ANOVA, ANCOVA, MANCOVA. In a section of my result, I did correlation of
>my
>data (2 categirical factors with 2 levels: Quantity & Quality; 2
>dependent
>var: Irid.area & Casa.PC1, and 1 co-var: SL). But as some traits (here
>Irid.area) are significantly influenced by the covariate (standard
>length,
>SL), I need to use the partial correlation. I know how to calculate it
>with
>JMP, but as I used R to analyse all of my data (first time in my life)
>for
>this manuscript, can anyone help me to find a solution for this
>problem? I
>got some libraries to calculate it (e.g. ppcor, ggm, etc.), but none of
>them fit to my required analysis (fitting covariate and subset group)
>in
>the model.
>
>Any help will be very much appreciated.
>
> > ### Datafrmame> data1 <- read.csv(file.choose(),header=TRUE)
>#Partial correlation.csv> data1 Quantity Quality SL Irid.area
> Casa.PC1
>1 High Low 16.38 10.31 1.711739555
>2 High High 15.95 16.52 0.013383537
>3 High High 15.69 12.74 2.228490878
>4 High Low 14.76 9.80 1.554975833
>5 High Low 14.63 12.95 1.823767970
>6 High High 14.32 14.21 3.152059841
>7 High High 14.95 12.57 2.069265040
>8 High Low 15.37 13.55 1.886027422
>9 High Low 14.73 14.18 1.127440602
>10 High High 16.08 15.98 1.435563307
>11 High High 15.78 16.76 2.433261686
>12 High Low 15.22 12.12 0.927454986
>13 High Low 14.22 10.91 2.328899576
>14 High High 14.47 11.03 1.522923487
>15 High Low 13.98 10.03 2.342535074
>16 High Low 14.99 11.44 0.749529924
>17 High High 16.51 20.16 2.993905677
>18 High High 14.83 16.82 2.227315597
>19 High Low 15.17 19.21 1.685063793
>20 High Low 16.29 20.31 1.551704440
>21 High High 16.23 15.03 1.982319336
>22 High High 14.18 14.80 1.839910851
>23 High Low 16.11 12.92 1.443240647
>24 High Low 13.95 7.60 2.034192171
>25 High High 17.54 17.80 2.188306237
>26 High Low 16.24 19.29 1.531264746
>27 High High 14.79 12.98 1.465644134
>28 High Low 15.87 14.85 1.372494892
>29 High High 16.09 13.71 1.462037152
>30 High Low 14.34 13.53 1.365588960
>31 High High 14.93 12.91 0.729212386
>32 High High 15.89 16.98 0.136175317
>33 High Low 16.11 11.93 1.442761666
>34 High Low 15.25 15.49 0.834442777
>35 High High 15.84 17.65 1.471713978
>36 High High 15.61 18.00 1.949457500
>37 High Low 15.42 13.87 0.200098471
>38 High Low 14.91 11.23 0.981988071
>39 High High 15.69 5.74 -0.445941360
>40 High High 15.13 9.07 1.387947896
>41 High Low 15.04 15.87 1.480980400
>42 High Low 17.08 17.24 2.620029423
>43 High High 15.85 12.47 0.027278890
>44 High High 15.35 10.44 2.597373230
>45 High Low 15.62 12.11 0.030653396
>46 High High 17.96 17.50 1.544922124
>47 High Low 17.25 17.87 1.705053951
>48 High Low 15.56 19.72 1.688867665
>49 High High 16.27 13.15 0.111371757
>50 High Low 16.68 15.43 1.538012366
>51 High High 15.78 15.07 0.744555741
>52 Low High 14.72 13.34 -0.682505420
>53 Low Low 14.93 14.07 -1.641494605
>54 Low High 13.94 13.22 -1.172268647
>55 Low High 14.01 18.65 -0.996656064
>56 Low Low 14.33 17.16 -1.789728167
>57 Low Low 14.57 12.43 -0.827526343
>58 Low High 14.01 15.29 -1.350691602
>59 Low Low 14.22 16.98 -1.688278221
>60 Low High 13.45 14.40 -1.182117327
>61 Low High 13.44 16.57 -1.358976542
>62 Low Low 14.76 15.58 0.334534454
>63 Low Low 14.85 17.65 0.251766383
>64 Low High 13.42 10.99 -0.526634460
>65 Low High 14.07 16.88 -1.112579922
>66 Low Low 14.15 16.41 -0.971918177
>67 Low Low 14.78 11.95 -1.179074800
>68 Low High 14.84 17.62 -0.777057705
>69 Low High 15.16 14.09 -1.224388816
>70 Low Low 14.60 15.03 -0.775478528
>71 Low High 13.74 10.01 -0.917153842
>72 Low High 13.54 12.34 -0.822895877
>73 Low Low 14.04 11.86 0.002789116
>74 Low High 15.73 18.50 -1.209469875
>75 Low Low 15.14 16.85 -0.479090055
>76 Low Low 14.86 17.32 -1.897204235
>77 Low High 14.43 11.20 0.469569392
>78 Low Low 14.01 15.55 -1.025059269
>79 Low High 14.20 11.67 -0.770451072
>80 Low High 16.16 17.34 -0.274527631
>81 Low Low 14.63 13.52 -1.070187945
>82 Low Low 15.83 14.85 -1.627211162
>83 Low High 14.70 14.81 -1.694118608
>84 Low High 13.91 14.48 -1.635459183
>85 Low Low 13.95 16.05 -1.449612666
>86 Low Low 14.03 12.58 -1.685968841
>87 Low High 14.82 13.57 -0.097426417
>88 Low High 14.32 12.16 -1.403512009
>89 Low Low 14.33 7.66 -1.336654713
>90 Low Low 15.01 10.15 -1.257019268
>91 Low High 14.01 9.79 -0.715404495
>92 Low Low 14.25 17.38 -1.296954022
>93 Low High 14.55 16.11 -0.616895943
>94 Low High 13.98 11.49 -0.654017365
>95 Low Low 15.59 8.43 -1.708330027
>96 Low Low 15.02 16.88 -1.352913634
>97 Low High 13.99 9.64 -0.499793618
>98 Low Low 13.98 12.25 -1.265336955
>99 Low High 13.94 13.79 -0.263925513
>100 Low Low 15.03 20.39 -0.720121308
>101 Low Low 13.93 14.63 -0.908570400> ### COrrelation
>test according to group> library("MASS")> with(data1,
cor.test(~
>Irid.area + Casa.PC1, subset=(Quantity=="High")))# gives cor,
df+2,
>p-values
> Pearson's product-moment correlation
>
>data: Irid.area and Casa.PC1
>t = 1.5795, df = 49, p-value = 0.1206
>alternative hypothesis: true correlation is not equal to 0
>95 percent confidence interval:
> -0.05905155 0.46734855
>sample estimates:
> cor
>0.2201142
>> with(data1, cor.test(~ Irid.area + Casa.PC1,
>subset=(Quantity=="Low")))# gives cor, df+2, p-values
> Pearson's product-moment correlation
>
>data: Irid.area and Casa.PC1
>t = -0.4275, df = 48, p-value = 0.6709
>alternative hypothesis: true correlation is not equal to 0
>95 percent confidence interval:
> -0.3342116 0.2205349
>sample estimates:
> cor
>-0.06159377
>> #### Effect size from two-way ANOVA ####> anova<- aov(Irid.area ~
>Quantity*Quality+SL, data=data1)> summary(anova) Df Sum
>Sq Mean Sq F value Pr(>F)
>Quantity 1 0.0 0.04 0.004 0.947
>Quality 1 0.3 0.26 0.032 0.859
>SL 1 149.5 149.49 18.027 5.03e-05 ***
>Quantity:Quality 1 0.2 0.18 0.022 0.883
>Residuals 96 796.1 8.29
>---
>Signif. codes: 0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1 >
>etaSquared( anova ) # effect size eta.sq
>eta.sq.part
>Quantity 5.682119e-02 0.0632542041
>Quality 6.577118e-05 0.0000781554
>SL 1.510031e-01 0.1521470868
>Quantity:Quality 1.922552e-04 0.0002284211> ### partial correlation
>(pcor) tests:> library(ggm)Loading required package: graphError in
>loadNamespace(i[[1L]], c(lib.loc, .libPaths())) :
> there is no package called ?BiocGenerics?In addition: Warning
>messages:1: package ?ggm? was built under R version 2.15.3 2: package
>?graph? was built under R version 3.0.1 Error: package ?graph? could
>not be loaded> data2<- data1[, c("Irid.area",
"Casa.PC1", "SL"),
>Quantity == "High"]> pcor(c("Irid.area",
"Casa.PC1",
>"SL"),var(data2))Error in match.arg(method) : 'arg' must
be NULL or a
>character vector> pc<-pcor(data2)> pc$estimate
> Irid.area Casa.PC1 SL
>Irid.area 1.0000000 -0.1313475 0.3387663
>Casa.PC1 -0.1313475 1.0000000 0.5061438
>SL 0.3387663 0.5061438 1.0000000
>
>$p.value
> Irid.area Casa.PC1 SL
>Irid.area 0.0000000000 1.896426e-01 3.647247e-04
>Casa.PC1 0.1896425857 0.000000e+00 6.258573e-09
>SL 0.0003647247 6.258573e-09 0.000000e+00
>
>$statistic
> Irid.area Casa.PC1 SL
>Irid.area 0.000000 -1.311637 3.564375
>Casa.PC1 -1.311637 0.000000 5.809698
>SL 3.564375 5.809698 0.000000
>
>$n
>[1] 101
>
>$gp
>[1] 1
>
>$method
>[1] "pearson"
>
>
>
>
>Cheers,
Dear All, Sorry to bother you again. As my previous mail was messy to understand, please find it again to give me a solution. I'd like to do a partial correaltion test ['pcor.test ()' or 'parcor()'] between Irid.area and Casa.PC1 variables controlling the influence of SL (co-variate) according to categorical group factors (Quantity or Quality group). Can anyone give me any example or source of package that will be appreciated. I already tried with "ggm" (doesn't work in R now), "ppcor", but failed to modify or adjust with my data to anlyse as I want. My effort:> data1 Quantity Quality SL Irid.area Casa.PC11 High Low 16.38 10.31 1.71173956 2 High High 15.95 16.52 0.01338354 3 High High 15.69 12.74 2.22849088 4 High Low 14.76 9.80 1.55497583 5 High Low 14.63 12.95 1.82376797 6 High High 14.32 14.21 3.15205984 7 High High 14.95 12.57 2.06926504 8 High Low 15.37 13.55 1.88602742 9 High Low 14.73 14.18 1.12744060 10 High High 16.08 15.98 1.43556331 11 Low Low 13.95 16.05 -1.44961267 12 Low Low 14.03 12.58 -1.68596884 13 Low High 14.82 13.57 -0.09742642 14 Low High 14.32 12.16 -1.40351201 15 Low Low 14.33 7.66 -1.33665471 16 Low Low 15.01 10.15 -1.25701927 17 Low High 14.01 9.79 -0.71540450 18 Low Low 14.25 17.38 -1.29695402 19 Low High 14.55 16.11 -0.61689594 20 Low High 13.98 11.49 -0.65401736> ### Correlation test according to group> library("MASS")> with(data1, cor.test(~ Irid.area + Casa.PC1, subset=(Quantity=="High")))# gives cor, df+2, p-valuesPearson's product-moment correlation data: Irid.area and Casa.PC1 t = -1.0507, df = 8, p-value = 0.3241 alternative hypothesis: true correlation is not equal to 0 95 percent confidence interval: -0.8020154 0.3604104 sample estimates: cor -0.3482398> with(data1, cor.test(~ Irid.area + Casa.PC1, subset=(Quantity=="Low")))# gives cor, df+2, p-valuesPearson's product-moment correlation data: Irid.area and Casa.PC1 t = 0.1209, df = 8, p-value = 0.9068 alternative hypothesis: true correlation is not equal to 0 95 percent confidence interval: -0.6031431 0.6547228 sample estimates: cor 0.04269795> ### ppcor tests:> y.data <- data.frame(+ v1=data1$Irid.area,+ v2=data1$Casa.PC1,+ f1=data1$Quantity,+ f2=data1$Quality)> library(ppcor)> pcor(y.data)Error in pcor(y.data) : 'x' must be numeric ## I also tried it with numeric values, but didn't give me the results same as JMP> # partial correlation between "v1" and "v2" given "f1" and "f2"> pcor.test(y.data$v1,y.data$v2,y.data[,c("f1","f2")])Error: is.numeric(y) || is.logical(y) is not TRUE???Then I tried with> #####pcor with ggm#####> library(ggm)# suggested by Andy Field et al. on their book "Discovering statistics using R".Loading required package: graphError: package ‘graph’ could not be loaded> examData2<- data1[, c("Irid.area", "Casa.PC1", "SL")]> maleExam<-subset(data1, Quantity == "High", select= c("Irid.area", "Casa.PC1"))> femaleExam<-subset(data1, Quantity == "Low", select= c("Irid.area", "Casa.PC1"))> cor(maleExam) Irid.area Casa.PC1Irid.area 1.0000000 -0.3482398 Casa.PC1 -0.3482398 1.0000000> cor(femaleExam) Irid.area Casa.PC1 Irid.area 1.00000000 0.04269795 Casa.PC1 0.04269795 1.00000000> # partial correlation between two var.> pcor.test(examData2$Irid.area, examData2$Casa.PC1, examData2$SL..COVARIATE.,Quantity == "High")Error in pcor.test(examData2$Irid.area, examData2$Casa.PC1, examData2$SL..COVARIATE., : 'use' should be either "rec" or "mat"!In addition: Warning messages:1: In if (use == "mat") { : the condition has length > 1 and only the first element will be used2: In if (use == "rec") { : the condition has length > 1 and only the first element will be used So, again can anyone help me to find out the solution, and sorry in advance to disturb you with the same issue. Cheers, Jewel e-mail:moshiurku@yahoo.com g-mail: mrahmankufmrt@gmail.com [[alternative HTML version deleted]]