# A thing
2

# Use print to show multiple things in the same cell
# Note that you can use single or double quotes for strings
print 2
print 'hello'

# Things can be stored as variables
a = 2
b = 'hello'
c = True  # This is case sensitive
print a, b, c

# The type function tells us the type of thing we have
print type(a)
print type(b)
print type(c)

# Standard math operators work as expected in numerical expressions
print 2 + 3
print 2 * 3
print 2 ** 3

# We can also assign numbers to named variables for use in expressions
a = 2
b = 3
print a + b
print a * b
print a ** b

# But be careful with integer division
# Also watch out for exponentiation - you want **, not ^
print 2 / 3
print 2 / 3.0

print 2 ** 3
print 2 ^ 3

# There are also operators for strings
print 'hello' + 'world'
print 'hello' * 3
#print 'hello' / 3

# Boolean operators compare two things
a = (1 > 3)
b = (3 == 3)
print a
print b
print a or b
print a and b

# There are thousands of functions that operate on things
print type(3)
print len('hello')
print round(3.3)

#round?
#round(
round(3.14159, 2)

# Many useful functions are in external packages
# Let's meet numpy
import numpy as np

# To see what's in a package, type the name, a period, then hit tab
#np?
#np.

# Some examples of numpy functions and "things"
print np.sqrt(4)
print np.pi  # Not a function, just a variable
print np.sin(np.pi)

# A string is actually an object
a = 'hello, world'
print type(a)

# Objects have bundled methods
#a.
print a.capitalize()
print a.replace('l', 'X')

r = 0.6
K = 100
n0 = 10

n1 = round(n0 + r*n0*(1 - n0/K))
n2 = round(n1 + r*n1*(1 - n1/K))

print n0, n1, n2
print type(n2)

print 'n2 more than 20: ', n2 > 20
print 'n2 is an integer: ', n2.is_integer()

# Lists are created with square bracket syntax
a = ['hi', 'hello', 'yo']
print a, type(a)

# Lists (and all collections) are also indexed with square brackets
# NOTE: The first index is zero, not one
print a[0]
print a[1]

# Lists can be sliced by putting a colon between indexes
# NOTE: The end value is not inclusive
print a[0:2]

# You can leave off the start or end if desired
print a[:2]
print a[2:]
print a[:]
print a[:-1]

# Lists are objects, like everything else, and have methods such as append
a.append('hiya')
print a

a.append([1,2])
print a

a.pop()
print a

# A very common use of lists is to create a range of integers
# This is easy in Python using the built-in range function like
# range(start, end); note that like the 'slice' notation the end value
# is *not* inclusive
print range(0, 10)

######################################
# This code deletes our old variables
try: del n0, n1, n2, r, K
except: pass
######################################

r = 0.6
K = 100
n = []
n.append(10)  # Append n0 in the first location

n.append(round(n[0] + r*n[0]*(1 - n[0]/K)))  # Append n1
n.append(round(n[1] + r*n[1]*(1 - n[1]/K)))  # Append n2
n.append(round(n[2] + r*n[2]*(1 - n[2]/K)))  # Append n3

print n

print "Grew by a factor of ", n[3]/n[0]  # or n[-1]/n[0]

# Make a dictionary of model parameters
params = {'n0': 10, 'r': 0.5}
print params
print params['r']

params['K'] = 200
print params

#print params['c']

K = 200
locals()['K']

# Make an array from a list
alist = [2, 3, 4]
blist = [5, 6, 7]
a = np.array(alist)
b = np.array(blist)
print a, type(a)
print b, type(b)

# Do arithmetic on arrays
print a**2
print np.sin(a)
print a * b
print a.dot(b), np.dot(a, b)

# Boolean operators work on arrays too, and they return boolean arrays
print a > 2
print b == 6

c = a > 2
print c
print type(c)
print c.dtype

# Indexing arrays
print a[0:2]

c = np.random.rand(3,3)
print c
print '\n'
print c[1:3,0:2]

c[0,:] = a
print '\n'
print c

# Arrays can also be indexed with other boolean arrays
print a
print b
print a > 2
print a[a > 2]
print b[a > 2]

b[a == 3] = 77
print b

# ndarrays have attributes in addition to methods
#c.
print c.shape
print c.prod()

# There are handy ways to make arrays full of ones and zeros
print np.zeros(5), '\n'
print np.ones(5), '\n'
print np.identity(5), '\n'

E = np.empty(5)
E.fill(2)

# You can also easily make arrays of number sequences
print np.arange(0, 10, 2)

r = 0.6
K = 100
n = np.zeros(100)
n[0] = 10

n[1] = round(n[0] + r*n[0]*(1 - n[0]/K))
n[2] = round(n[1] + r*n[1]*(1 - n[1]/K))
n[3] = round(n[2] + r*n[2]*(1 - n[2]/K))

print n

t = np.arange(0, 25, 0.25)
t = np.linspace(0, 25, 100, endpoint=False)
print t

print n[t == 0.5]

# A basic for loop - don't forget the white space!
wordlist = ['hi', 'hello', 'bye']
for word in wordlist:
    print word + '!'

# Indentation error: Fix it!
for word in wordlist:
    new_word = word.capitalize()
   print new_word + '!' # Bad indent

# Sum all of the values in a collection using a for loop
numlist = [1, 4, 77, 3]

total = 0
for num in numlist:
    total = total + num
    
print "Sum is", total

# Often we want to loop over the indexes of a collection, not just the items
print wordlist

wordrange = np.arange(len(wordlist))
print wordrange

for i in wordrange:
    print i, wordlist[i]

# While loops are useful when you don't know how many steps you will need,
# and want to stop once a certain condition is met.
step = 0
prod = 1
while prod < 100:
    step = step + 1
    prod = prod * 2
    print step, prod
    
print 'Reached a product of', prod, 'at step number', step

# Load up pylab, a useful plotting library
%matplotlib inline
from matplotlib import pyplot as plt

# Make some x and y data and plot it
y = np.arange(100)**2
plt.plot(y)

r = 0.6
K = 100
n = np.zeros(100, dtype=float)
n[0] = 10

steps = range(1, 100)
for i in steps:
    n[i] = round(n[i-1] + r*n[i-1]*(1 - n[i-1]/K))
    
plt.plot(n)

r = 0.6
K = 100
n = np.zeros(100)
n[0] = 10

i = 1
while n[i-1] < 90:
    n[i] = round(n[i-1] + r*n[i-1]*(1 - n[i-1]/K))
    print i, n[i]
    i = i + 1
    
print 'Ended with a population size of', n[i-1], 'at step', i-1
print 'Population sizes to this step were', n[0:i]

# A simple if statement
x = 3
if x > 0:
    print 'x is positive'
elif x < 0:
    print 'x is negative'
else:
    print 'x is zero'

# If statements can rely on boolean variables
x = -1
test = (x > 0)
print type(test); print test

if test:
    print 'Test was true'

# Set up initial model parameters
r = 0.6
K = 100
n = np.zeros(100)
n[0] = 10

# Loop through all time steps
steps = np.arange(1, 100)
for i in steps:
    cata = (np.random.rand() < 0.1)  # Random catastrophe 10% of time
    if cata:
        n[i] = n[0]
    else:
        n[i] = round(n[i-1] + r*n[i-1]*(1 - n[i-1]/K))
    
plt.plot(n)

print 'Population was above 50 in', np.sum(n > 50), 'out of 100 steps'

# We've been using functions all day
x = 3.333333
print round(x, 2)
print np.sin(x)

# It's very easy to write your own functions
def multiply(x, y):
    return x*y

# Once a function is "run" and saved in memory, it's available just like any other function
print type(multiply)
print multiply(4, 3)

# It's useful to include docstrings to describe what your function does
def say_hello(time, people):
    '''
    Function says a greeting. Useful for engendering goodwill
    '''
    return 'Good ' + time + ', ' + people

say_hello('afternoon', 'friends')

# All arguments must be present, or the function will return an error
#say_hello('afternoon')

# Keyword arguments can be used to make some arguments optional by giving them a default value
# All mandatory arguments must come first, in order
def say_hello(time, people='friends'):
    return 'Good ' + time + ', ' + people

say_hello('afternoon')

say_hello('afternoon', 'students')

def logistic_growth(r, K, n0, p=0.1):
    '''
    Function to simulate discrete time stochastic logistic growth

    Arguments
    ---------
    r : float
        Reproductive rate of population
    K : float
        Carrying capacity
    n0 : float
        Initial population
    p : float
        Probability of a catastrophe resetting the population to n0

    Returns
    -------
    n : ndarray
        Array of 100 time steps of population values
    '''
    
    # Set up population array
    n = np.zeros(100, dtype=float)
    n[0] = n0

    # Loop through all time steps
    steps = np.arange(1, 100)
    for i in steps:
        cat = (np.random.rand() < p)  # Random catastrophe 10% of time
        if cat:
            n[i] = 10
        else:
            n[i] = grow_one_step(n[i-1], r, K)
    
    return n

def grow_one_step(nold, r, K):
    '''Calculate new population from old one.'''
    return round(nold + r*nold*(1 - nold/K))

a1 = logistic_growth(0.6, 100, 10)
a2 = logistic_growth(0.6, 100, 10)
a3 = logistic_growth(0.6, 100, 10)

plt.hold(True)
plt.plot(a1)
plt.plot(a2)
plt.plot(a3)

# import pop
# reload(pop) # This line forces the re-import of pop if you make changes and run this cell again
# plot(pop.logistic_growth(0.6, 100, 10))