<jupyter_start><jupyter_text>Homework 10: Linear Regression **Reading**: * [Linear Regression](https://www.inferentialthinking.com/chapters/15/2/Regression_Line.html)* [Method of Least Squares](https://www.inferentialthinking.com/chapters/15/3/Method_of_Least_Squares.html)* [Least Squares Regression](https://www.inferentialthinking.com/chapters/15/4/Least_Squares_Regression.html) Please complete this notebook by filling in the cells provided. Directly sharing answers is not okay, but discussing problems with the course staff or with other students is encouraged. Refer to the policies page to learn more about how to learn cooperatively.For all problems that you must write our explanations and sentences for, you **must** provide your answer in the designated space. Moreover, throughout this homework and all future ones, please be sure to not re-assign variables throughout the notebook! For example, if you use `max_temperature` in your answer to one question, do not reassign it later on.<jupyter_code># Don't change this cell; just run it. 

import numpy as np
from datascience import *

# These lines do some fancy plotting magic.
import matplotlib
%matplotlib inline
import matplotlib.pyplot as plt
plt.style.use('fivethirtyeight')
import warnings
warnings.simplefilter('ignore', FutureWarning)

import otter
grader = otter.Notebook()<jupyter_output><empty_output><jupyter_text>Exploring the PTEN Gene with Linear Regression  1. PTEN Linear Regression This week's homework is about linear regression. The dataset we'll be using is from the Cancer Cell Line Encyclopedia -- you can read more about this database in this [paper](https://www.nature.com/articles/s41586-019-1186-3) and interact with the data yourself at the online portal [here](https://portals.broadinstitute.org/ccle).The specific dataset we'll be taking a look at is expression data for the PTEN gene in around 1000 cell lines. The PTEN gene is a tumor-suppressing gene, and mutations in the PTEN gene are associated with many types of cancer. A cell line is group of cells that are kept alive and replicate indefinitely in culture (grown in petri dishes, for example).Run the following cell to load the `pten` table. The `pten` table has four columns, a column for the specific `Cell Line`, a column for the `Copy Number`, which is how many times a copy of a portion of the PTEN gene is found in the DNA of that cell line, `mRNA Expression (Affy)`, and `mRNA Expression (RNAseq)`.<jupyter_code># Just run this cell
pten = Table().read_table("pten.csv")
pten.show(5)
# Just run this cell
pten.hist("Copy Number", bins = np.arange(-1, 1.5, 0.5))<jupyter_output><empty_output><jupyter_text>**Question 1**Looking at the histogram above, we want to check whether or not `Copy Number` is in standard units. For this question, compute the mean and the standard deviation of the values in `Copy Number` and assign these values to `copy_number_mean` and `copy_number_sd` respectively. After you calculate these values, assign `is_su` to either `True` if you think that `Copy Numbers` is in standard units or `False` if you think otherwise.<!--BEGIN QUESTIONname: q1_1manual: false--><jupyter_code>copy_number = pten.column("Copy Number")
copy_number_mean = ...
copy_number_sd = ...
is_su = ...
print(f"Mean: {copy_number_mean}, SD: {copy_number_sd}, Is in standard units?: {is_su}")
grader.check("q1_1")<jupyter_output><empty_output><jupyter_text>**Question 2**Create the function `standard_units` so that it converts the values in the array `arr` to standard units. We'll then use `standard_units` to create a new table, `pten_su`, that converts all the values in the table `pten` to standard units.<!--BEGIN QUESTIONname: q1_2manual: false--><jupyter_code>def standard_units(arr):
    ...

# DON'T DELETE OR MODIFY ANY OF THE LINES OF CODE BELOW IN THIS CELL
pten_su = Table().with_columns("Cell Line", pten.column("Cell Line"),
                               "Copy Number SU", standard_units(pten.column("Copy Number")),
                               "mRNA Expression (Affy) SU", standard_units(pten.column("mRNA Expression (Affy)")),
                               "mRNA Expression (RNAseq) SU", standard_units(pten.column("mRNA Expression (RNAseq)"))                             
                              )
pten_su.show(5)
grader.check("q1_2")<jupyter_output><empty_output><jupyter_text>You should always visually inspect your data, before numerically analyzing any relationships in your dataset. Run the following cell in order to look at the relationship between the variables in our dataset.<jupyter_code># Just run this cell
pten_su.scatter("Copy Number SU", "mRNA Expression (Affy) SU")
pten_su.scatter("Copy Number SU", "mRNA Expression (RNAseq) SU")
pten_su.scatter("mRNA Expression (Affy) SU", "mRNA Expression (RNAseq) SU")<jupyter_output><empty_output><jupyter_text>**Question 3**Which of the following relationships do you think has the highest correlation (i.e. highest absolute value of `r`)? Assign `highest_correlation` to the number corresponding to the relationship you think has the highest correlation.1. Copy Number vs. mRNA Expression (Affy)2. Copy Number vs. mRNA Expression (RNAseq)3. mRNA Expression (Affy) vs. mRNA Expression (RNAseq)<!--BEGIN QUESTIONname: q1_3manual: false--><jupyter_code>highest_correlation = ...
grader.check("q1_3")<jupyter_output><empty_output><jupyter_text>**Question 4**Now, using the `standard units` function, define the function `correlation` which computes the correlation between `arr1` and `arr2`.<!--BEGIN QUESTIONname: q1_4manual: false--><jupyter_code>def correlation(arr1, arr2):
    ...

# This computes the correlation between the different variables in pten
copy_affy = correlation(pten.column("Copy Number"), pten.column("mRNA Expression (Affy)"))
copy_rnaseq = correlation(pten.column("Copy Number"), pten.column("mRNA Expression (RNAseq)"))
affy_rnaseq = correlation(pten.column("mRNA Expression (Affy)"), pten.column("mRNA Expression (RNAseq)"))

print(f" \
      Copy Number vs. mRNA Expression (Affy) Correlation: {copy_affy}, \n \
      Copy Number vs. mRNA Expression (RNAseq) Correlation: {copy_rnaseq}, \n \
      mRNA Expression (Affy) vs. mRNA Expression (RNAseq) Correlation: {affy_rnaseq}")
grader.check("q1_4")<jupyter_output><empty_output><jupyter_text>**Question 5**If we switch what we input as arguments to `correlation`, i.e. found the correlation between `mRNA Expression (Affy)` vs. `Copy Number` instead of the other way around, would the correlation change? Assign `correlation_change` to either `True` if you think yes, or `False` if you think no.<!--BEGIN QUESTIONname: q1_5manual: false--><jupyter_code>correlation_change = ...
grader.check("q1_5")<jupyter_output><empty_output><jupyter_text>**Question 6**Looking at both the scatter plots after Question 2 and the correlations computed in Question 4, describe a pattern you see in the relationships between the variables.<!--BEGIN QUESTIONname: q1_6manual: true--> _Type your answer here, replacing this text._ **Question 7**Let's look at the relationship between mRNA Expression (Affy) vs. mRNA Expression (RNAseq) only. Define a function called `regression_parameters` that returns the parameters of the regression line as a two-item array containing the slope and intercept of the regression line as the first and second elements respectively. The function `regression_line` takes in two arguments, an array of `x` values, and an array of `y` values.<!--BEGIN QUESTIONname: q1_7manual: false--><jupyter_code>def regression_parameters(x, y):
    ...
    slope = ...
    intercept = ...
    return make_array(slope, intercept)

parameters = regression_parameters(pten.column("mRNA Expression (Affy)"), pten.column("mRNA Expression (RNAseq)"))
parameters
grader.check("q1_7")<jupyter_output><empty_output><jupyter_text>**Question 8**If we switch what we input as arguments to `regression_parameters`, i.e. found the parameters for the regression line for `mRNA Expression (RNAseq)` vs. `mRNA Expression (Affy)` instead of the other way around, would the regression parameters change (would the slope and/or intercept change)? Assign `parameters_change` to either `True` if you think yes, or `False` if you think no.<!--BEGIN QUESTIONname: q1_8manual: false--><jupyter_code>parameters_change = ...
grader.check("q1_8")<jupyter_output><empty_output><jupyter_text>**Question 9**Now, let's look at how the regression parameters look like in standard units. Use the table `pten_su` and the function `regression_parameters`, and assign `parameters_su` to a two-item array containing the slope and the intercept of the regression line for mRNA Expression (Affy) in standard units vs. mRNA Expression (RNAseq) in standard units.<!--BEGIN QUESTIONname: q1_9manual: false--><jupyter_code>parameters_su = ...
parameters_su
grader.check("q1_9")<jupyter_output><empty_output><jupyter_text>**Question 10**Looking at the array `parameters_su`, what do you notice about the slope and intercept values specifically? Relate them to another value we already calculated in a previous question, as well as relate them to an equation.<!--BEGIN QUESTIONname: q1_10manual: true--> _Type your answer here, replacing this text._ **Question 11**The oldest and most commonly used cell line in Biology is the HeLa cell line, named after <NAME>, whose cervical cancer cells were taken without her consent in 1951 to create this cell line. The issue of data privacy and consent is very important to data science! You can read more about this topic [here](https://www.hopkinsmedicine.org/henriettalacks/).The HeLa cell line is missing from our dataset. If we know that the HeLa mRNA Expression (Affy) value is 8.2, what is the predicted mRNA Expression (RNAseq) value? Use the values in `parameters` that we derived in Question 1.7, and assign the result to `hela_rnaseq`.<!--BEGIN QUESTIONname: q1_11manual: false--><jupyter_code>hela_rnaseq = ...
hela_rnaseq
grader.check("q1_11")<jupyter_output><empty_output><jupyter_text>**Question 12**Compute the predicted mRNA Expression (RNAseq) values from the mRNA Expression (Affy) values in the `pten` table. Use the values in the `parameters` array from Question 1.7, and assign the result to `predicted_rnaseq`. We'll plot your computed regression line with the scatter plot from after question 1.2 of mRNA Expression (Affy) vs. mRNA Expression (RNAseq).<!--BEGIN QUESTIONname: q1_12manual: true--><jupyter_code>predicted_rnaseq = ...

# DON'T CHANGE/DELETE ANY OF THE BELOW CODE IN THIS CELL
(pten.with_column("Predicted mRNA Expression (RNAseq)", predicted_rnaseq)
 .select("mRNA Expression (Affy)", "mRNA Expression (RNAseq)", "Predicted mRNA Expression (RNAseq)")
 .scatter("mRNA Expression (Affy)"))
plt.ylabel("mRNA Expression (RNAseq)");<jupyter_output><empty_output><jupyter_text>Fitting a least-squares regression line Recall that the least-square regression line is the unique straight line that minimizes root mean squared error (RMSE) among all possible fit lines. Using this property, we can find the equation of the regression line by finding the pair of slope and intercept values that minimize root mean squared error. **Question 13**Define a function called `RMSE`. It should take two arguments:1. the slope of a line (a number)2. the intercept of a line (a number).It should return a number that is the root mean squared error (RMSE) for a line defined with the arguments slope and intercept used to predict mRNA Expression (RNAseq) values from mRNA Expression (Affy) values for each row in the `pten` table.*Hint: Errors are defined as the difference between the actual `y` values and the predicted `y` values.**Note: if you need a refresher on RMSE, here's the [link](https://www.inferentialthinking.com/chapters/15/3/Method_of_Least_Squares.htmlRoot-Mean-Squared-Error) from the textbook*<!--BEGIN QUESTIONname: q1_13manual: false--><jupyter_code>def RMSE(slope, intercept):
    affy = pten.column("mRNA Expression (Affy)")
    rnaseq = pten.column("mRNA Expression (RNAseq)")
    predicted_rnaseq = ...
    ...

# DON'T CHANGE THE FOLLOWING LINES BELOW IN THIS CELL
rmse_example = RMSE(0.5, 6)
rmse_example
grader.check("q1_13")<jupyter_output><empty_output><jupyter_text>**Question 14**What is the RMSE of a line with slope 0 and intercept of the mean of `y` equal to?*Hint 1: The line with slope 0 and intercept of mean of `y` is just a straight horizontal line at the mean of `y`**Hint 2: What does the formula for RMSE become if we input our predicted `y` values in the formula. Try writing it out on paper! It should be a familiar formula.*<!--BEGIN QUESTIONname: q1_14manual: true--> _Type your answer here, replacing this text._ **Question 15**Find the parameters that minimizes RMSE of the regression line for mRNA Expression (Affy) vs. mRNA Expression (RNAseq). Assign the result to `minimized_parameters`.If you haven't tried to use the `minimize` [function](http://data8.org/sp20/python-reference.html) yet, now is a great time to practice. Here's an [example from the textbook](https://www.inferentialthinking.com/chapters/15/3/Method_of_Least_Squares.htmlnumerical-optimization).*Hint: Use the `RMSE` function in Question 1.13***NOTE: When you use the minimize function, please pass in `smooth=True` as the second argument to this function. You'll need to do this, otherwise, your answer will be incorrect**<!--BEGIN QUESTIONname: q1_15manual: false--><jupyter_code>minimized_parameters = ...
minimized_parameters
grader.check("q1_15")<jupyter_output><empty_output><jupyter_text>**Question 16**The slope and intercept pair you found in Question 1.15 should be very similar to the values that you found in Question 1.7. Why were we able to minimize RMSE to find the same slope and intercept from the previous formulas?<!--BEGIN QUESTIONname: q1_16manual: true--> _Type your answer here, replacing this text._ **Question 17**If we had instead minimized mean squared error (MSE), would we have gotten the same slope and intercept of the minimized root mean squared error (RMSE) results? Assign `same_parameters` to either `True` if you think yes, or `False` if you think no.<!--BEGIN QUESTIONname: q1_17manual: false--><jupyter_code>same_parameters = ...
same_parameters
grader.check("q1_17")<jupyter_output><empty_output><jupyter_text>Let's look at the scatter plot of the relationship between mRNA Expression (Affy) and mRNA Expression (RNAseq) again:<jupyter_code>pten.scatter("mRNA Expression (Affy)", "mRNA Expression (RNAseq)")<jupyter_output><empty_output><jupyter_text>**Question 18**Using a linear regression model, would we be able to obtain accurate predictions for most of the points? Explain why or why not.<!--BEGIN QUESTIONname: q1_18manual: true--> _Type your answer here, replacing this text._  2. Properties of Binary Distributions Binary distributions arise in regular everyday life, and as data scientists you will encounter them constantly. A binary distribution is a distribution across two categories: such as voting in support of a proposition or voting against it on your local ballot, flipping heads or tails, having heart disease or not having heart disease. Generally we represent 'yes' or `True` as 1, and 'no' or `False` as 0. Binary distributions have some special properties that make working with them especially easy! The intent of this section of the homework is to walk you through these properties, so we decided to make all of the tests for this section public (i.e. there are no hidden tests to worry about for this section only).  Question 1Let's generate a random binary distribution of 0's and 1's. Assign `binary_options` to the correct array of possible values in a binary distribution (i.e. look at the previous sentence).<!--BEGIN QUESTIONname: q2_1manual: false--><jupyter_code>binary_options = ...

# DON'T DELETE/MODIFY ANY OF THE CODE IN THIS CELL BELOW
sample_size = 100
binary_sample = np.random.choice(binary_options, sample_size)

# Run this to see a histogram of this random distribution.
Table().with_columns("Value", make_array(1, 0), "Number in Sample", make_array(sum(binary_sample), sample_size - sum(binary_sample))).barh("Value")
grader.check("q2_1")<jupyter_output><empty_output><jupyter_text>Question 2The first property you should note is that the proportion of ones in a binary distribution is equal to the mean of the distribution. [Think about why this is true](https://www.inferentialthinking.com/chapters/14/1/Properties_of_the_Mean.htmlProportions-are-Means). Complete the following cell to show that this is the case for your `binary_sample`. Assign `number_of_ones` and `number_of_zeros` to the number of `1`'s and the number of `0`'s respectively from your `binary_sample`.<!--BEGIN QUESTIONname: q2_2manual: false--><jupyter_code>number_of_ones = ...
number_of_zeros = ...


# DON'T DELETE/MODIFY ANY OF THE CODE BELOW IN THIS CELL
number_values = len(binary_sample)
sum_of_binary_sample = sum(binary_sample)
# Remember that the mean is equal to the sum divided by the number of items
mean_binary_sample = sum_of_binary_sample / number_values

# Don't change this!
print(f"In your binary sample there were {number_of_ones} ones and {number_of_zeros} zeros. 1*{number_of_ones} + 0*{number_of_zeros} = {number_of_ones}")
print(f"The sum of values in your sample was {sum_of_binary_sample}, divided by the number of items, {number_values}, gives us a mean of {mean_binary_sample}")
print(f"The proportion of ones in your sample was {number_of_ones} ones, divided by the number of items, {number_values}, gives us a value of {mean_binary_sample}" )
print('Those values are equal!')
grader.check("q2_2")<jupyter_output><empty_output><jupyter_text>Since the proportion of ones is the same as the mean, the Central Limit Theorem applies! That is, if we resample our sample a lot of times, the distribution of the proportion of ones in our resamples will be roughly normal, with a predictable center and spread!<jupyter_code># Just run this cell
resampled_proportion_of_ones = make_array()

for i in np.arange(5000):
    resample = Table().with_column("Value", binary_sample).sample()
    resample_proportion_ones = resample.where("Value", 1).num_rows / resample.num_rows
    resampled_proportion_of_ones = np.append(resampled_proportion_of_ones, resample_proportion_ones)
    
Table().with_column('Resampled Proportions', resampled_proportion_of_ones).hist()<jupyter_output><empty_output><jupyter_text>Let's generate a table where each row has a different number of ones and zeros that we'll use for the following parts.<jupyter_code># Just run this cell
possible_number_ones = np.arange(sample_size + 1)
possible_number_zeros = sample_size - possible_number_ones

possibilities_table = Table().with_columns("Values of One", possible_number_ones, "Values of Zero", possible_number_zeros)
possibilities_table.show(5)<jupyter_output><empty_output><jupyter_text>Question 3The second important property of binary distributions is that the standard deviation of every binary distribution is equal to:$$\sqrt{\text{proportion_ones} *\text{proportion_zeros}}$$While this property is useful in some cases, a more useful extension of this property is that it tells us that the maximum standard deviation for a binary distribution is 0.5!Let's explore why that is the case!Complete the `binary_std_formula` function below so that it returns the standard deviation of a binary distribution according to the formula above.<!--BEGIN QUESTIONname: q2_3manual: false--><jupyter_code>def binary_std_formula(row):
    num_ones = row.item("Values of One")
    num_zeros = row.item("Values of Zero")
    
    sum_ones_and_zeros = ...
    prop_ones = ...
    prop_zeros = ...
    ...

# DON'T DELETE/MODIFY ANY OF THE LINES BELOW IN THIS CELL
possibilities_table = possibilities_table.with_column("Formula SD", possibilities_table.apply(binary_std_formula))
possibilities_table.show(5)
grader.check("q2_3")<jupyter_output><empty_output><jupyter_text>Here's another function that takes in a row object from a table, generates a sample that has the same number of ones and zeros as the row specifies, and then returns the standard deviation of that table. You should be able to understand exactly what this function does! It also does the same thing as above, where we return the standard deviation, but we just use `np.std` for this function.<jupyter_code># Just run this cell
def binary_std(row):
    values = make_array()
    for i in np.arange(row.item("Values of One")):
        values = np.append(values, 1)
    for i in np.arange(row.item("Values of Zero")):
        values = np.append(values, 0)
    return np.std(values)

possibilities_table = possibilities_table.with_column("Empirical SD", possibilities_table.apply(binary_std))
possibilities_table.show(5)<jupyter_output><empty_output><jupyter_text>All the values are the same! Let's see what this formula means!<jupyter_code># Just run this cell
possibilities_table.scatter("Values of One", "Formula SD")<jupyter_output><empty_output><jupyter_text>What a beautiful curve!Looking at that curve, we can see that maximum value is $0.5$, which occurs in the middle of the distribution, when the two categories have equal proportions (proportion of ones = proportion of zeros = $\frac{1}{2}$).  (OPTIONAL, NOT IN SCOPE) Logarithmic Plots A kind of visualization you will frequently encounter as a data scientist is a scatter plot or line plot that uses a logarithmic scale. This **Optional** section will cover how to read and generate logarithmic plots. Since this is optional, there is no autograded/free response questions for these sections. Just read, run cells, and explore. What is a logarithm? A logarithm helps us find the inverse of an equation that uses exponentials. Specifically, if$$a^y = x$$Then$$\log_a{x} = y$$The most commonly used $a$, which is known as the base of the logarithm, is $e$, which is equivalent to about 2.718, or 10 (for powers of 10).We can use `numpy` to take logs in Python! By default, np.log uses a base of e.<jupyter_code>make_array(np.log(np.e), np.log(np.e**2), np.log(100))<jupyter_output><empty_output><jupyter_text>Back to the visualization: when we are plotting trends that grow exponentially, such as the line$$ y = e^x$$our y-axis needs to have a large range of values, which makes it difficult to understand.Let's see what this looks like:<jupyter_code>x = np.arange(0, 10, 1/100)
y = 10 ** x

Table().with_columns("X", x, "Y", y).scatter(0,1)<jupyter_output><empty_output><jupyter_text>Note that since $10^{10}$ is so big, we can't really see what's happening at all to the y values when they have x values below 8.One solution to this to change our y and/or x axis so that instead of having even spaces between the tick marks, our marks grow by an uneven factor. We do this by making the tick marks go on a logarithmic scale, and we'll then be able to understand our data better!<jupyter_code>Table().with_columns("X", x, "Y", y).scatter(0,1)
plt.yscale("log")<jupyter_output><empty_output><jupyter_text>Now we can tell what's happening to the y values for every x value!Note how the y values start at $10^0=1$, and increase by a *factor* of $10$ each mark - the next mark is $10^1 = 10$, then $10^2=100$.You still read this plot like a normal plot, so at a value of $x=5, y=10^5=10000$.How do you calculate intermediate values? At a value like $x = 2.5$ it looks like the y value is somewhere in-between $10^1$ and $10^3$. In this graph with a logarithmic scale, you would say that $y=10^{2.5} \approx 316$. When visualizing data about the spread of diseases, you will commonly run into plots with logarithmic scales, such as this example from the New York Times. Make sure to always know what the scales of the data are! Image is from https://www.nytimes.com/2020/03/20/health/coronavirus-data-logarithm-chart.html ---To double-check your work, the cell below will rerun all of the autograder tests.<jupyter_code>grader.check_all()<jupyter_output><empty_output><jupyter_text>SubmissionMake sure you have run all cells in your notebook in order before running the cell below, so that all images/graphs appear in the output. The cell below will generate a zip file for you to submit. **Please save before exporting!**<jupyter_code># Save your notebook first, then run this cell to export your submission.
grader.export(pdf=False)<jupyter_output><empty_output>
StarcoderdataJupyterStructuredCleanDedup