Correlation and Dependency¶

Finding links between variables is a very important task. The main purpose of data science is to find relationships between variables, and to understand how these relationships can help us make better decisions.

Machine learning models are also sensitive to the number of variables and how they relate and affect each other, so finding correlations and dependencies can help us make better use of our machine learning algorithms.

Let’s use the Telco Churn dataset to understand how we can find links between different variables in VerticaPy.

import verticapy as vp

churn = vp.read_csv("customers.csv")
churn.head(100)

	Abc customerID Varchar(20) 100%	...	123 TotalCharges Numeric(9,3) 99%	010 Churn Boolean 100%
1	0002-ORFBO	...	593.3	❌
2	0003-MKNFE	...	542.4	❌
3	0004-TLHLJ	...	280.85	✅
4	0017-DINOC	...	2460.55	❌
5	0020-JDNXP	...	1993.2	❌
6	0022-TCJCI	...	2791.5	✅
7	0023-XUOPT	...	1215.6	✅
8	0040-HALCW	...	1090.6	❌
9	0042-JVWOJ	...	471.85	❌
10	0064-SUDOG	...	224.5	❌
11	0080-EMYVY	...	727.85	❌
12	0083-PIVIK	...	5567.55	❌
13	0089-IIQKO	...	3767.4	❌
14	0093-EXYQL	...	3673.6	❌
15	0093-XWZFY	...	4036.85	✅
16	0103-CSITQ	...	6252.7	❌
17	0111-KLBQG	...	2861.45	❌
18	0115-TFERT	...	2317.1	✅
19	0117-LFRMW	...	1448.8	✅
20	0133-BMFZO	...	181.65	✅
21	0141-YEAYS	...	2401.05	❌
22	0156-FVPTA	...	1152.7	✅
23	0164-XAIRP	...	470.2	❌
24	0178-CIIKR	...	58.0	❌
25	0187-WZNAB	...	1972.35	❌
26	0188-GWFLE	...	33.7	❌
27	0206-OYVOC	...	864.2	❌
28	0218-QNVAS	...	7113.75	❌
29	0219-QAERP	...	576.65	❌
30	0221-NAUXK	...	219.5	❌
31	0233-FTHAV	...	4765.0	❌
32	0237-YFUTL	...	5405.8	❌
33	0254-FNMCI	...	7624.2	❌
34	0254-KCJGT	...	4354.45	❌
35	0254-WWRKD	...	3431.75	❌
36	0259-GBZSH	...	181.5	✅
37	0265-PSUAE	...	1522.7	❌
38	0266-GMEAO	...	8058.55	❌
39	0277-ORXQS	...	3364.55	❌
40	0282-NVSJS	...	355.9	❌
41	0307-BCOPK	...	326.65	❌
42	0310-MVLET	...	6010.05	✅
43	0320-JDNQG	...	2331.3	✅
44	0322-YINQP	...	48.55	✅
45	0324-BRPCJ	...	6851.65	✅
46	0330-IVZHA	...	330.15	✅
47	0345-XMMUG	...	4854.3	❌
48	0348-SDKOL	...	5885.4	✅
49	0354-WYROK	...	2911.3	✅
50	0357-NVCRI	...	471.7	❌
51	0363-QJVFX	...	3432.9	✅
52	0365-BZUWY	...	1742.5	❌
53	0369-ZGOVK	...	1992.2	❌
54	0376-OIWME	...	3366.05	❌
55	0376-YMCJC	...	1943.2	✅
56	0378-NHQXU	...	1460.65	✅
57	0379-DJQHR	...	5398.6	❌
58	0380-NEAVX	...	3343.15	❌
59	0383-CLDDA	...	5897.4	❌
60	0386-CWRGM	...	967.85	❌
61	0388-EOPEX	...	139.4	✅
62	0396-YCHWO	...	3474.2	❌
63	0411-EZJZE	...	170.5	❌
64	0412-UCCNP	...	3175.85	❌
65	0415-MOSGF	...	44.4	✅
66	0420-HLGXF	...	4036.0	❌
67	0422-OHQHQ	...	295.95	❌
68	0423-UDIJQ	...	447.9	❌
69	0431-APWVY	...	2598.95	✅
70	0436-TWFFZ	...	5714.2	❌
71	0439-IFYUN	...	1294.6	❌
72	0461-CVKMU	...	1900.25	✅
73	0464-WJTKO	...	1460.85	❌
74	0479-HMSWA	...	2715.3	❌
75	0480-BIXDE	...	1743.05	❌
76	0481-SUMCB	...	4735.35	❌
77	0484-JPBRU	...	996.45	❌
78	0506-LVNGN	...	349.65	✅
79	0513-RBGPE	...	2278.75	❌
80	0533-UCAAU	...	4140.1	❌
81	0537-QYZZN	...	1857.75	❌
82	0547-HURJB	...	235.8	❌
83	0549-CYCQN	...	3974.7	❌
84	0565-JUPYD	...	6590.8	❌
85	0567-GGCAC	...	438.9	❌
86	0575-CUQOV	...	5867.0	❌
87	0576-WNXXC	...	2510.2	✅
88	0577-WHMEV	...	1374.9	❌
89	0581-MDMPW	...	2072.75	❌
90	0601-WZHJF	...	667.7	✅
91	0602-DDUML	...	3894.4	❌
92	0603-TPMIB	...	1534.05	❌
93	0620-XEFWH	...	84.2	❌
94	0621-HJWXJ	...	5029.05	❌
95	0628-CNQRM	...	1544.05	✅
96	0635-WKOLD	...	2921.75	❌
97	0637-YLETY	...	1555.65	✅
98	0650-BWOZN	...	1359.45	❌
99	0655-YDGFJ	...	1323.7	❌
100	0657-DOGUM	...	2985.25	❌

The Pearson correlation coefficient is a very common correlation function. In this case, it helped us to find linear links between the variables. Having a strong Pearson relationship means that the two input variables are linearly correlated.

churn.corr(method = "pearson")

We can see that tenure is well-correlated to the TotalCharges, which makes sense.

churn.scatter(["tenure", "TotalCharges"])

churn.corr(["tenure", "TotalCharges"], method = "pearson")
Out[1]: 0.825880460933202

Note, however, that having a low Pearson relationship imply that the variables aren’t correlated. For example, let’s compute the Pearson correlation coefficient between tenure and TotalCharges to the power of 20.

churn["TotalCharges^20"] = churn["TotalCharges"] ** 20

churn.corr(["tenure", "TotalCharges^20"], method = "pearson")
Out[3]: 0.224994408804537

We know that the tenure and TotalCharges are strongly linearly correlated. However we can notice that the correlation between the tenure and TotalCharges to the power of 20 is not very high. Indeed, the Pearson correlation coefficient is not robust for monotonic relationships, but rank-based correlations are. Knowing this, we’ll calculate the Spearman’s rank correlation coefficient instead.

churn.corr(method = "spearman", show = False)

	...	"Churn"	"TotalCharges^20"
"SeniorCitizen"	...	0.150889328176473	0.105795342303725
"Partner"	...	-0.150447544959177	0.343930553215626
"Dependents"	...	-0.164221401579725	0.0866797760484616
"tenure"	...	-0.369620778763435	0.883103368818293
"PhoneService"	...	0.0119419800290031	0.0838048547856037
"PaperlessBilling"	...	0.191825331666468	0.151669712799097
"MonthlyCharges"	...	0.184839285783758	0.633958405301206
"TotalCharges"	...	-0.233211018585104	1.0
"Churn"	...	1.0	-0.233211018585104
"TotalCharges^20"	...	-0.233211018585104	1.0

churn.corr(method = "spearman")

The Spearman’s rank correlation coefficient determines the monotonic relationships between the variables.

churn.corr(["tenure", "TotalCharges^20"], method = "spearman")
Out[4]: 0.883103368818293

We can notice that Spearman’s rank correlation coefficient stays the same if one of the variables can be expressed using a monotonic function on the other. The same applies to Kendall rank correlation coefficient.

churn.corr(method = "kendall")

Notice that the Kendall rank correlation coefficient will also detect the monotonic relationship.

churn.corr(["tenure", "TotalCharges^20"], method = "kendall")
Out[5]: 0.731699318287362

However, the Kendall rank correlation coefficient is very computationally expensive, so we’ll generally use Pearson and Spearman when dealing with correlations between numerical variables.

Binary features are considered numerical, but this isn’t technically accurate. Since binary variables can only take two values, calculating correlations between a binary and numerical variable can lead to misleading results. To account for this, we’ll want to use the Biserial Point method to calculate the Point-Biserial correlation coefficient. This powerful method will help us understand the link between a binary variable and a numerical variable.

churn.corr(method = "biserial")

Lastly, we’ll look at the relationship between categorical columns. In this case, the Cramer's V method is very efficient. Since there is no position in the Euclidean space for those variables, the Cramer's V coefficients cannot be negative (which is a sign of an opposite relationship) and they will range in the interval [0,1].

churn.corr(method = "cramer")

Sometimes, we just need to look at the correlation between a response and other variables. The parameter focus will isolate and show us the specified correlation vector.

churn.corr(method = "cramer", focus = "Churn")

Sometimes a correlation coefficient can lead to incorrect assumptions, so we should always look at the coefficient p-value.

churn.corr_pvalue("Churn", "customerID", method = "cramer")
Out[6]: (nan, nan)

We can see that churning correlates to the type of contract (monthly, yearly, etc.) which makes sense: you would expect that different types of contracts differ in flexibility for the customer, and particularly restrictive contracts may make churning more likely.

The type of internet service also seems to correlate with churning. Let’s split the different categories to binaries to understand which services can influence the global churning rate.

churn["InternetService"].one_hot_encode()
churn.corr(
    method = "spearman",
    focus = "Churn",
    columns = [
        "InternetService_DSL",
        "InternetService_Fiber_optic",
    ],
)

We can see that the Fiber Optic option in particular seems to be directly linked to a customer’s likelihood to churn. Let’s compute some aggregations to find a causal relationship.

churn["contract"].one_hot_encode()
churn.groupby(
    [
        "InternetService_Fiber_optic",
    ],
    [
        "AVG(tenure) AS tenure",
        "AVG(totalcharges) AS totalcharges",
        'AVG("contract_month-to-month") AS "contract_month-to-month"',
        'AVG("monthlycharges") AS "monthlycharges"',
    ],
)

	123 InternetService_Fiber_optic Integer 100%	...	123 contract_month-to-month Float(22) 100%	123 monthlycharges Float(22) 100%
1	0	...	0.442614644033443	43.7882442361287
2	1	...	0.68733850129199	91.5001291989664

It seems that users with the Fiber Optic option tend more to churn not because of the option itself, but probably because of the type of contracts and the monthly charges the users are paying to get it. Be careful when dealing with identifying correlations! Remember: correlation doesn’t imply causation!

Another important type of correlation is the autocorrelation. Let’s use the amazon dataset to understand it.

from verticapy.datasets import load_amazon

amazon = load_amazon()
amazon.head(100)

	📅 date Date 100%	...	Abc state Varchar(32) 100%	123 number Int 100%
1	1998-01-01	...	AMAPÁ	0
2	1998-01-01	...	AMAZONAS	0
3	1998-01-01	...	DISTRITO FEDERAL	0
4	1998-01-01	...	ESPÍRITO SANTO	0
5	1998-01-01	...	MARANHÃO	0
6	1998-01-01	...	PARANÁ	0
7	1998-01-01	...	PIAUÍ	0
8	1998-01-01	...	RORAIMA	0
9	1998-01-01	...	SERGIPE	0
10	1998-01-01	...	SÃO PAULO	0
11	1998-02-01	...	GOIÁS	0
12	1998-02-01	...	MATO GROSSO DO SUL	0
13	1998-02-01	...	MINAS GERAIS	0
14	1998-02-01	...	PARAÍBA	0
15	1998-02-01	...	SANTA CATARINA	0
16	1998-02-01	...	SÃO PAULO	0
17	1998-03-01	...	AMAPÁ	0
18	1998-03-01	...	BAHIA	0
19	1998-03-01	...	MATO GROSSO DO SUL	0
20	1998-03-01	...	PARÁ	0
21	1998-03-01	...	PERNAMBUCO	0
22	1998-03-01	...	RIO GRANDE DO SUL	0
23	1998-04-01	...	CEARÁ	0
24	1998-04-01	...	PARANÁ	0
25	1998-04-01	...	PARAÍBA	0
26	1998-04-01	...	PARÁ	0
27	1998-05-01	...	ALAGOAS	0
28	1998-05-01	...	AMAPÁ	0
29	1998-05-01	...	MARANHÃO	0
30	1998-05-01	...	MATO GROSSO DO SUL	0
31	1998-05-01	...	RIO GRANDE DO SUL	0
32	1998-05-01	...	TOCANTINS	0
33	1998-06-01	...	ESPÍRITO SANTO	6
34	1998-06-01	...	RIO DE JANEIRO	3
35	1998-06-01	...	RIO GRANDE DO NORTE	1
36	1998-06-01	...	SÃO PAULO	451
37	1998-07-01	...	ESPÍRITO SANTO	37
38	1998-07-01	...	MARANHÃO	274
39	1998-07-01	...	MATO GROSSO	360
40	1998-07-01	...	MATO GROSSO DO SUL	3712
41	1998-07-01	...	PARAÍBA	0
42	1998-07-01	...	PARÁ	638
43	1998-07-01	...	RONDÔNIA	365
44	1998-07-01	...	SÃO PAULO	596
45	1998-08-01	...	ALAGOAS	1
46	1998-08-01	...	BAHIA	815
47	1998-08-01	...	DISTRITO FEDERAL	48
48	1998-08-01	...	ESPÍRITO SANTO	38
49	1998-08-01	...	MARANHÃO	1176
50	1998-08-01	...	MATO GROSSO	228
51	1998-08-01	...	MINAS GERAIS	875
52	1998-08-01	...	PIAUÍ	711
53	1998-08-01	...	RIO GRANDE DO SUL	9
54	1998-08-01	...	RORAIMA	0
55	1998-08-01	...	SERGIPE	0
56	1998-09-01	...	AMAPÁ	20
57	1998-09-01	...	DISTRITO FEDERAL	33
58	1998-09-01	...	PIAUÍ	1991
59	1998-09-01	...	RORAIMA	2
60	1998-10-01	...	AMAZONAS	83
61	1998-10-01	...	GOIÁS	1034
62	1998-10-01	...	MATO GROSSO	576
63	1998-10-01	...	PARAÍBA	179
64	1998-10-01	...	PARÁ	3665
65	1998-10-01	...	PIAUÍ	2586
66	1998-10-01	...	SERGIPE	0
67	1998-11-01	...	ALAGOAS	19
68	1998-11-01	...	AMAPÁ	131
69	1998-11-01	...	CEARÁ	575
70	1998-11-01	...	DISTRITO FEDERAL	0
71	1998-11-01	...	MARANHÃO	2237
72	1998-11-01	...	RIO DE JANEIRO	6
73	1998-11-01	...	RIO GRANDE DO SUL	28
74	1998-11-01	...	SÃO PAULO	488
75	1998-12-01	...	BAHIA	82
76	1998-12-01	...	MARANHÃO	1399
77	1998-12-01	...	MATO GROSSO	100
78	1998-12-01	...	PARAÍBA	51
79	1998-12-01	...	PERNAMBUCO	59
80	1998-12-01	...	RIO DE JANEIRO	1
81	1998-12-01	...	RONDÔNIA	33
82	1998-12-01	...	TOCANTINS	9
83	1999-01-01	...	ALAGOAS	58
84	1999-01-01	...	GOIÁS	14
85	1999-01-01	...	MATO GROSSO	239
86	1999-01-01	...	MINAS GERAIS	36
87	1999-01-01	...	PARÁ	87
88	1999-01-01	...	PERNAMBUCO	102
89	1999-01-01	...	RONDÔNIA	1
90	1999-01-01	...	SÃO PAULO	7
91	1999-01-01	...	TOCANTINS	36
92	1999-02-01	...	ACRE	0
93	1999-02-01	...	CEARÁ	16
94	1999-02-01	...	MATO GROSSO	69
95	1999-02-01	...	MATO GROSSO DO SUL	28
96	1999-02-01	...	PERNAMBUCO	13
97	1999-02-01	...	RONDÔNIA	1
98	1999-02-01	...	SANTA CATARINA	2
99	1999-02-01	...	TOCANTINS	1
100	1999-03-01	...	AMAPÁ	2

Our goal is to predict the number of forest fires in Brazil. To do this, we can draw an autocorrelation plot and a partial autocorrelation plot.

amazon.acf(
    column = "number",
    ts = "date",
    by = ["state"],
    p = 24,
    method = "pearson",
)

amazon.pacf(
    column = "number",
    ts = "date",
    by = ["state"],
    p = 8,
)

We can see the seasonality forest fires.

It’s mathematically impossible to build the perfect correlation function, but we still have several powerful functions at our disposal for finding relationships in all kinds of datasets.