<jupyter_start><jupyter_text>Numerical v. NumericalWe've already seen how to explore a single numerical variable and it makes sense that each numerical variable should be explored on its own. What we're interested in now is how to look at two different variables and see if there might be a relationship between them.We've already seen that we can look at two descriptive statistics that measure the relationships between two numeric variables. Pearson's $r$ (often just called "correlationc coefficient") is a measure of the positive or negative *linear* relationship between two variables. Spearman's $\rho$ is the measure of the increasing or decreasing monotonic relationship between each variable.<jupyter_code>random.seed(57483249)

data = pd.DataFrame(np.concatenate([stats.norm.rvs( 90.0, 9.0, 100), 
                                    stats.norm.rvs(120.0, 12.0, 100)]))
data.columns = ["X"]
data["Y"] = data["X"].apply(lambda x: 5.0 + 0.4 * x + 
                            stats.norm.rvs( 0.0, 8.0, 1)[0])
data.head()<jupyter_output><empty_output><jupyter_text>Although we can assume that each variable will have been examined individually, let's go ahead and calculate the descriptive statistics for each variable:<jupyter_code>data.describe()<jupyter_output><empty_output><jupyter_text>Depending on how we arrived at this point, we may or may not already know the correlation between the two variables. We want to look at both the Pearson correlation coefficient and the Spearman correlation coefficient. As we remarked earlier, Pearson measures *linear* relationships whereas Spearman measures the *monotonicity* of the relationship.It's also worth noting that the Coefficient of Variation ($R^2$) is just this Pearson correlation coefficient, squared.<jupyter_code>print( "r   =", stats.pearsonr(data["X"], data["Y"])[0])
print( "rho =", stats.spearmanr(data["X"], data["Y"])[0])<jupyter_output>r   = 0.6989968396163851
rho = 0.7046786169654242<jupyter_text>There is a fairly strong correlation between the two variables. The convention is that $r$ must be greater than 0.5 for there to be an "interesting" correlation. However, don't forget that $r$ only measures the degree of the linear relationship. Scatter PlotsLet's visualize the data using a scatter plot, the typical visualization for two numerical variables:<jupyter_code>figure = plt.figure(figsize=(10, 6))

axes = figure.add_subplot(1, 1, 1)
axes.scatter( data[ "X"], data["Y"], marker="o", color="DimGray")
axes.set_title("X v. Y")
axes.set_xlabel( "X")
axes.set_ylabel( "Y")

plt.show()
plt.close()<jupyter_output><empty_output><jupyter_text>When visualizing data it's important to be cognizant of two things: measurement artifacts and overstriking. Measurement artifacts arise when some variable is measured on a fixed scale or convention, for example, miles per gallon is almost always reported in integral values. What would it look like if $Y$ was measured on an integral scale?<jupyter_code>figure = plt.figure(figsize=(10, 6))

integral_y = data["Y"].apply(lambda x: np.floor(x/5) * 5)

axes = figure.add_subplot(1, 1, 1)
axes.scatter( data[ "X"], integral_y, marker="o", color="DimGray")
axes.set_title("X v. Y")
axes.set_xlabel( "X")
axes.set_ylabel( "Y")

plt.show()
plt.close()<jupyter_output><empty_output><jupyter_text>For example, suppose $Y$ was measured in increments of 5 (5, 10, 15, 20). As we can see above, this creates a striating pattern that can dominate ability to find larger patterns in the data.There is also overstriking when many observations have the same or similar values.The solution to the first problem is to add a bit of "jitter" to the data to get rid of the striping effect--but not so much that we end up with false values. The solution to the second problem is to apply an alpha to the markers so that if there is overstriking, those points will be darker as more points are added.<jupyter_code>figure = plt.figure(figsize=(10, 6))

integral_y = data["Y"].apply(lambda x: np.floor(x/5) * 5)
jittered_y = integral_y + stats.norm.rvs( 0, 0.5, integral_y.size)

axes = figure.add_subplot(1, 1, 1)
axes.scatter( data["X"], jittered_y, \
             marker="o", color="DimGray", alpha=0.7)

axes.set_title("X v. Y")
axes.set_xlabel("X")
axes.set_ylabel("Y")
plt.show()
plt.close()<jupyter_output><empty_output><jupyter_text>Even with just some transparency added to the markers, we can start to see some kind of density arise in the plot. You can think of this as a third dimension to the data. The jitter also breaks up the striations a bit; this probably isn't enough jitter (it's normally distributed "noise" with mean = 0 and standard deviation 0.5). Because the distances between the bins is 5 units, we could probably use a jitter of 1-1.5. This would keep most points from "jumping" into another bin. I left the jitter small so you could see it being applied. HexbinningThere may be cases where you have so much data that neither jitter nor transparency will help. In that case, you can use hexbinning. Hexbinning divides the XY-plane into hexagons of a specified size that are then used a bins in a 2d histogram (viewed from the top). A picture is worth a 1,000 words (in both senses).Let's make some "big" data:<jupyter_code>big_data = pd.DataFrame(np.concatenate([stats.norm.rvs( 90.0, 9.0, 10000), \
                                        stats.norm.rvs(120.0, 12.0, 10000)]))
big_data.columns = ["X"]
big_data["Y"] = big_data["X"].apply(lambda x: 5.0 + 0.4 * x \
                                    + stats.norm.rvs( 0.0, 8.0, 1)[0])
big_data.head()<jupyter_output><empty_output><jupyter_text>And plot it the usual way:<jupyter_code>figure = plt.figure(figsize=(10, 6))

axes = figure.add_subplot(1, 1, 1)
axes.scatter(big_data[ "X"], big_data["Y"], marker="o", color="DimGray")
axes.set_title("X v. Y")
axes.set_xlabel( "X")
axes.set_ylabel( "Y")

plt.show()
plt.close()<jupyter_output><empty_output><jupyter_text>Yikes. That's not helpful. Let's try hexbinning:<jupyter_code>from matplotlib import cm

figure = plt.figure(figsize=(10, 6))

axes = figure.add_subplot(1, 1, 1)
axes.hexbin(big_data["X"], big_data["Y"], gridsize=20, cmap=cm.Greys)

plt.show()
plt.close()<jupyter_output><empty_output><jupyter_text>The hexbinning reveals a lot more structure to the data. And just like with other types of histograms, you may need to try different "grid" sizes. Another approach would be to *sample* your data.  TrendlinesWe may also want to add a trendline to our data based on the outcome of the Pearson or Spearman coefficient. At this stage we generally go for a non-parametric approach like LOESS (locally weighted regression). We'll have more to say about regression in later chapters.`matplotlib` doesn't have Lowess out of the box, but we can use `statsmodels`:<jupyter_code>#https://stackoverflow.com/questions/36252434/predicting-on-new-data-using-locally-weighted-regression-loess-lowess/37060980

import statsmodels.api as sm

lowess = sm.nonparametric.lowess(data["Y"], data["X"], frac=.3)
lowess_x = list(zip(*lowess))[0]
lowess_y = list(zip(*lowess))[1]

figure = plt.figure(figsize=(10, 6))

axes = figure.add_subplot(1, 1, 1)
axes.scatter( data[ "X"], data["Y"], marker="o", color="DimGray")
axes.plot(lowess_x, lowess_y, color="DarkRed")

axes.set_title("Plot of X v. Y with LOWESS")
axes.set_xlabel( "X")
axes.set_ylabel( "Y")

plt.show()
plt.close()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.11013868318298102, 'fraction_numerical': 0.02598626632556887, 'mean_word_length': 4.844217151848937, 'pattern_counts': {'":': 0, '<': 40, '<?xml version=': 0, '>': 40, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '10', 'dir': 'jupyter-structured-clean-dedup', 'id': '661629', 'license': 'MIT', 'n_tokens_mistral': 2546, 'n_tokens_neox': 2292, 'n_words': 873, 'path': 'fundamentals_2018.9/eda/matrix/numerical_v_numerical.ipynb', 'repo_name': 'topseer/APL_Great', 'size': '7427'}