مهملترین کاری که ممکنه این هفته با پایتون، پردازش تصویر، آر و آمار استنباطی انجام بدین: بررسی کبریت‌های توکلی

حکایتی هست در این باره که راکفلکر شروع کارش از یک کارخونه کبریت سازی بوده. این داستان مدعی است که این مولتی میلیاردر که زمان جوانی در کارخونه کبریت سازی کار می کرده، روزی پیش مدیرش می ره و می گه ایده ای داره برای چند برابر کردن درآمد کارخونه و در مقابل این قول که در بخشی از این درآمد سهیم بده ایده اش رو می ده:

جعبه کبریت یک کشویی داره که کبریت ها اون تو هستن. این رو برعکس توی جعبه بذارین. هر کس جعبه رو باز می کنه بیشتر کبریت ها به زمین می ریزن و طرف خیلی زودتر از حالتی که کشویی درست نصب شده باشه، مجبور می شه یک جعبه کبریت دیگه بخره. چون کبریت چیز ارزونی است کسی به این امر توجه نخواهد کرد و هیچ کس هم شاکی نخواهد شد و شما چند برابر بیشتر جعبه کبریت خواهید فروخت.

ایده خبیثانه جالبیه ولی چرا اصولا از اول کبریت های کمتری توی جعبه نذاریم؟ اصولا اگر می شه اینجوری پول درآورد مطمئنا کشور ما جای بسیار مناسبی براشه.

من به کبریت مشکوکم

برای بررسی این مساله به مغازه می ریم و یک بسته ده تایی کبریت توکلی که ظاهرا تنها بازیگر بزرگ کبریت در ایرانه رو می خریم. اونو به خونه می‌یاریم. ملحفه سبزمون رو رو زمین پهن می کنیم و به عنوان یک تنبل واقعی، فقط قسمتی که قراره تو عکس باشه رو اتو می کنیم:

بسته ده تایی کبریت

و بعد می ریم سراغ ده تا جعبه کبریتی که خریدیم و پشتشون مدعی است که تعداد متوسط کبریت‌های هر جعبه چهل تا است:

کبریت توکلی

هر جعبه کبریت رو جدا جدا روی قسمت اتو شده و صاف خالی کرده، عکس می گیریم:

تصویر اولیه کبریت ها

اینکار رو برای هر ده تا جعبه کبریت انجام می دیم و اتاق پر از کبریت رو به قصد نشستن پشت کامپیوتر ترک میکنیم.

جعبه های کبریت

اول نیازمند یک برنامه هستیم که بتونه در یک عکس کبریت ها رو بشمره. بعد با یک اسکریپت این برنامه رو برای همه عکس ها اجرا خواهیم کرد. من برای اینکار اول کتابخونه‌های cv2 رو امتحان کردم ولی با numpy نتیجه بهتر و قشنگ تری گرفتم. در قدم اول کافیه عکس رو بخونیم، اون ۱ به کتابخونه می گه عکس رو بعد از خوندن خاکستری کنه که کار ما رو راحتتر می کنه:

kebrit = scipy.misc.imread(fileName, 1) # gray-scale image

و خروجی چیز شبیه این است:

تصویر کبریت قبل از پردازش تصویر

اینکه تصویر رنگی است به خاطر اینه که به شکل پیش فرض از پالت jet برای نمایش استفاده کرده ایم که تصویر رو به صورت حرارتی نشون می ده. می گن این برای کارهای علمی راحتتر از خاکستری است. به راحتی می بینین که چیزی که می خونیم با اینکه در سبز معمولی یا خاکستری به نظر یکدست می‌یاد ولی پر از تیزی و نقطه نقطه است. مثل هر عکس دوربینی دیگه. برای حل این مشکل تصویر رو اسموت می کنیم. اسموت کردن در پردازش تصویر معمولا با استفاده از فیلتر گاوس به دست می یاد که کانولوشن گرفتنی است بین این ورودی و تابع گاوس. البته نگران نباشین، من چون سه بار معادلات دیفرانسیل رو افتادم و آخرش هم یکی دیگه به جام امتحان داد، پایه ام قوی شده (:

kebrit_smooth = ndimage.gaussian_filter(kebrit, 6)

حالا متغیر kebrit_smooth یک تصویر نرم شده از تصویر اصلی است:

تصویر با فیلتر گاوسی پردازش تصویر

قدم بعدی اینه که من بیام و یک ترش‌هولد روی ورودی بذارم و بگم بیخیال هر چیزی بشه که از این عدد خاص کوچکتر است. این عدد با کمی تجربه و کمی سعی خطا به دست اومده. دستور طبیعی باید این می بود:

tresh = 120
labeled, objectsNum = ndimage.label(kebrit_smooth < tresh)

تابع ndimage.label ورودی خودش رو بررسی می کنه و تعداد اجسام به هم پیوسته توی اون رو می شمره. اما مشکل اینه که حتی تک پیکسل های منفرد هم شمرده می شن. برای جلوگیری از این امر من بهش یک ساختار می دم که شکلم حداقل لازمه چنین ساختاری داشته باشه:

tresh = 120
removeOnes = np.ones((3,3), dtype="bool8")
labeled, objectsNum = ndimage.label(kebrit_smooth < tresh, structure=removeOnes)

دقت کنین که اون استراکچر در اصل می گه حداقل تصویر قابل تشخیص من باید این باشه:

array([[ True,  True,  True],
       [ True,  True,  True],
       [ True,  True,  True]], dtype=bool)

و نتیجه اش این می شه که نقطه های منفرد دیگه شمرده نشن. البته انتظارم اینه که با اسموت کردن مرحله قبل این نقطه ها به حداقل رسیده باشن. نتیجه این ترکیب چنین تصویری است:

پردازش تصویر نهایی نوک کبریت ها

و خب چنین خروجی متنی ای:

10boxes/IMG_2608.JPG , 29

حالا کافیه یک اسکریپت ساده، همه ده تا فایل رو به اون برنامه بده:

#!/bin/bash

for f in *JPG
do
    ../count.py $f
done

و به خروجی زیر برسیم:

$ ./doall.sh 
IMG_2606.JPG , 35
IMG_2607.JPG , 30
IMG_2608.JPG , 29
IMG_2609.JPG , 32
IMG_2610.JPG , 29
IMG_2611.JPG , 33
IMG_2612.JPG , 36
IMG_2613.JPG , 38
IMG_2614.JPG , 37
IMG_2615.JPG , 38

وقتشه سراغ زبون مورد علاقمون R بریم. یک زبان تخصصی برای کارهای آماری و وررفتن با اعداد و ماتریس ها و رفیق رفقاشون. فایل رو می خونیم و بخش مورد نظر رو جدا می کنیم:

> tavakoli <- read.csv(file="results.csv",head=FALSE,sep=",")
> matches <- tavakoli[2]$V2

حالا اطلاعات هر جعبه کبریت رو دارم. بذارین یک نگاه سریع بهش بندازیم:

> summary (matches)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  29.00   30.50   34.00   33.70   36.75   38.00 

بعله. می بینیم که میزان متوسط پایینتر از ۳۴ است. یعنی تعداد متوسط این ده قوطی کبریت ۶ تا کمتر از چیزی است که ادعا شده یا به بیان صحیح‌تر ۱۵٪ کمتر به ما جنس فروختن. جالبیش اینه حتی یک جعبه از این ده جعبه هم به چیزی که ادعا می شد میانگین است نرسیده. بذارین جلوتر بریم:

> hist(matches, xlim=c(28,42), ylim=c(0,4), main="هیستوگرام تعداد کبریت های توکلی در هر جعبه", sub="www.jadi.net", breaks=8)
> curve(dnorm(x, mean=mean(matches), sd=sd(matches))*10, add=TRUE, col="red", lwd=2) 

و روی این نمودار با کمی آمار احتمالات می تونیم احتمال اینکه در یک جعبه کبریت اتفاقی توکلی که خریده‌ایم و مدعی داشتن ۴۰ کبریت است، چهل یا بیشتر کبریت وجود داشته باشد را حساب کنیم:

> 1-pnorm(40, mean=mean(matches), sd=sd(matches)) #یک منهای سمت چپ نمودار نرمال در نقطه ۴۰ کبریت
[1] 0.03970968

بله. با اینکه کبریت توکلی مدعی است در هر جعبه اش تقریبا ۴۰ کبریت وجود داره، برنامه پردازش تصویر و نمودارهای نرمال ما نشون می دن که احتمال اینکه واقعا در یک جعبه کبریت توکلی چهل یا بیشتر کبریت باشه، سه صدم درصد بیشتر نیست

آمار شیرین است

نتیجه‌ها

۰- موسسه استاندارد و حقوق مصرف کننده خاصی نداریم یا هنوز این مساله رو ندیدن
۱- راکفلر باید بیاد پیش تولید کننده‌های ما لنگ بندازه
۲- آمار شیرین و فان است
۳- می تونیم با کارهای علمی بامزه هم تفریح کنیم هم چیز یاد بگیریم

مرتبط

و البته مثل همیشه هر جاییش ممکنه اشتباه داشته باشه. خوشحال می شم دوستان حرفه ای اصلاح کنن یا توسعه بدن و به این فکر کنیم که کاش آمار رو اینطوری به ما درس می دادن (:

کارگاه‌های آموزشی پایتون و لینوکس گروه کاربران لینوکس کرج

اگر این هفته به برنامه دانشگاه امیر کبیر نمی‌رسین، می تونین هفته بعد سراغ کارگاه‌های آموزشی پایتون و لینوکس گروه کاربران کرج هم برین. روز بیست و پنجم کارگاه لینوکس است و روز بیست و ششم کارگاه پایتون. آزاده و رایگان و خوشحال.

poster_karaj

زمان کارگاه لینوکس: ۲۵ اردیبهشت ۱۳۹۳ ساعت ۱۰
زمان کارگاه پایتون: ۲۶ اردیبهشت ۱۳۹۳ ساعت ۱۰
مکان: کرج خیابان شهید بهشتی میدان حصارک دانشگاه خوارزمی

آپدیت: بر اساس قوانین دانشگاه ، برای حضور حتما باید ثبت نام کرده باشین… لطفا از این دو لینک پیش ثبت نام رو انجام بدین: http://events.karajlug.org/events/7 و http://events.karajlug.org/events/8

دعوت به کنفرانس پایتون

پنج سالی هست که عده ای از طرفداران پایتون دارن روی برگزاری یک کنفرانس کار می کنن. این کنفرانس اولین کنفرانس پایتون ایرانه و قراره در تاریخ ۱۰ آبان در فرهنگسرای مادر دانشگاه علم و صنعت توی تهران برگزار بشه. اگر علاقمند به شرکت هستین سری به سایت بزنین و اگر علاقمندتر هستین، تا اول‌ آبان مقاله‌تون رو بفرستین. منم با علاقه دوست دارم توش شرکت کنم و خوشبختانه جمعه هم هست.

چرا گنو لینوکس رو دوست دارم: دانلود کردن خودکار هر فصل جدید به محض انتشار

سعید پرسیده:

سلام جادی جان.
یه سوال داشتم ازت.
من چند مدته دارم آموزش سی شارپ رو یاد میگیرم از طریق سایت webtarget.ir
و اینجوریه که هر از چند مدت یه آموزش میزاره و من دنبال میکنم. میخ,استم بدونم راهی هست که من بدون اینکه مراجعه کنم به سایت و هر دفعه فایل های پی دی اف رو دانلود کنم بتونم با یه اسکریپت یا یه کد توی ترمینال تمام آموزهای سی شارپ رو توی یه فولدر دانلود کنم؟
ضمنا اینم آدرس یکی از فایل های پی دی اف آموزش سی شارپ هست و جاهایی که هر دفعه تغییر میکنه رو توی لینک دوم با ستاره برات مشخص میکنم.
http://dl.webtarget.ir/027-cSharpTime/session-33/cSharpTimeSession-33.pdf
http://dl.webtarget.ir/027-cSharpTime/session-**/cSharpTimeSession-**.pdf

خب.. سوال های اینجوری مثل یک پازل یا جدول بامزه هستن که جلوی تلویزیون لم می دیم و حلشون می کنیم. البته از نظر عقلی این سایت باید حتما RSS داشته باشه که بشه به سادگی دنبالش کرد ولی خب… ایده شما چیه؟ ایده من اینه که کافیه یکی از این فایل ها توی یک دایرکتوری باشه و ما روزی یکبار از اون دایرکتوری ls بگیریم و در آخرین فایل دانلود شده (بالاترین عدد در اسم ها) نگاه کنیم و عدد فایل رو جدا کنیم و بهش یکی اضافه کنیم و دانلودش کنیم. اگر فایل جدیدی بود دانلود می شه و اگر نبود چیزی دانلود نمی شه. فردا هم روز از نو روزی از نو.

برنامه ساده می شه این:

#!/usr/bin/python

import os, re

# یک ال.اس. می گیرم از فایل های مشابه اون کتاب و سورت و خط آخر رو جدا می کنم
f = os.popen('ls -1  cSharp*pdf | sort | tail -1')
lastdl = f.read()
# متغیر حاوی بخشی از اسم فایل است که عدد سریال توش قرار داره + ۱
lastNum = str(int(re.search ('cSharpTimeSession-(\d+).pdf', lastdl).group(1))+1)

# یک کامند می سازم که با دبلیوگِت شماره بعدی رو دانلود می کنه	
dlCommand = "wget http://dl.webtarget.ir/027-cSharpTime/session-" \
			+ lastNum + "/cSharpTimeSession-" \
			+ lastNum + ".pdf"
# اجراش می کنم
f = os.popen(dlCommand)

و البته اگر بخوایم به حالت های خاص جواب بدیم و مثلا خودمون اگر هیچ فایلی نبود به اولین فایل یک عدد بدیم (ظاهرا در سایتشون اولین نسخه شماره ۰۴ است) یا اگر زیر ۹ بودیم یک صفر اولش اضافه کنیم و … برنامه می شه این:

#!/usr/bin/python

import os, re

try:
	f = os.popen('ls -1  cSharp*pdf | sort | tail -1')
	lastdl = f.read()
	lastNum = int(re.search ('cSharpTimeSession-(\d+).pdf', lastdl).group(1))
except:
	lastNum = 3 #first file there is 4. so we'll assume that the current one is 3
	
if lastNum < 9:
	lastNum = "0" + str(lastNum+1)
else:
	lastNum = str (lastNum + 1)
	
dlCommand = "wget http://dl.webtarget.ir/027-cSharpTime/session-" \
			+ lastNum + "/cSharpTimeSession-" \
			+ lastNum + ".pdf"
f = os.popen(dlCommand)

منطقا این رو باید در یک کرون بذاریم یا سری اول با تکنیک watch python ./autodlcSharpbook.py دانلودش کنم که همه شماره ها رو بگیره یا مثلا اگر خروجی دستور موفقیت آمیز بود یک ایمیل بزنیم به صاحب جریان که فایل بهش اتچ باشه یا چنین چیزهایی. اونش با شما و سلیقه‌های شخصی‌تون.

پ.ن. اینم من حین نوشتن این برنامه. اصلا دلیلی که الان نوشتم این بود که دو تا مانیتور داشتم. فردا باید مانیتورم رو تحویل بدم ولی باید درخواست یک مانیتور اضافی بکنم. به نظرم هیچ چیز به اندازه دو - یا بیشتر - مانیتور کامپیوتر رو برای یک گیک لذت بخش نمی کنه. احتمالا عددی به اسم w وجود داره که اینطوری تعریف می شه «تعداد مانیتوری که یک سیستم ویندوزی باید داشته باشه تا یک علاقمند لینوکس با لذت پشتش بشینه و کار کنه»

jadi_while_programming

چرا گنو/لینوکس رو دوست دارم: شمردن کلمات به کار رفته شده در وبلاگ

در ادامه بررسی کلمات وبلاگ (بعد از شمردن استفاده از مهمل و کارهای آماری فوق العاده جالب روی اون) یک ایده دیگه که در کامنت‌ها مطرح شده رو پی‌می‌گیریم: شمردن کلمات و کاربرد اونها.

برنامه پایتون رو اینطوری تغییر می‌دم:

#!/usr/bin/python
# -*- coding: utf8 -*- 

from xml.dom import minidom
import xml.etree.cElementTree as et
import re

tree=et.parse('wordpress.2012-12-05.xml')
root=tree.getroot();

wordCount = {}

for child in root.iter('item'):
	date = child.find('wppost_date').text[:7] #find the year and month
	body = child.find('content_encoded').text # post content
	title = child.find('title').text 	  # post title

	try:
		fulltext = title + "\n" + body # all the text in the post = title + body
	except:
		pass

	fulltext = re.sub(ur'[_»«"\'&?؟a-zA-Z‌0-9/=.*+\n-%<>:;،؛,\-)(،۱۲۳۴۵۶۷۸۹۰]', ' ', fulltext) #replace extra chars
	words = fulltext.split() # words is a list of all words in this post

	for word in words:
		wordCount[word] = wordCount.get(word, 0) + 1 # wordCount[word]++

for word in wordCount:
	thisLine = (str(wordCount[word]) + "\t" + word)
	print thisLine.encode("utf-8", "ignore")


کد کاملا واضحه: هر پست رو نگاه می کنه، حروف غیرفارسی تیتر و متن رو با فاصله جایگزین می کنه و بعد تعداد کلمات رو جمع می‌زنه و همین روند رو روی تمام پست‌ها ادامه می‌ده. به عبارت دیگه خروجی چیزی شبیه به این خواهد بود:

...
4	همسرش
1	آکر
3	خرمش
1	ویسمن
2	خرما
1	یکباری
1	مانغو
2	احساسم
1	عصبي
61	رشد
1	رشت
1	تریلیان
5	همسرم
32	هیات
1	پورتال
2	پیشانی
6	مدیربسته
4	لری
4	وجدان
...


و البته مشخصه که خیلی طولانی‌تر. بذارین بشمریم که من کلا در زندگی چند کلمه استفاده کردم تو وبلاگم:

jadi@jeducation:~/Downloads/weblog_word_usage$ ./count_all_words.py | wc -l
27880


هوم.. بیست و هفت هزار و هشتصد و هشتاد کلمه (: بدک نیست (: چیز خوبیه که روی وبلاگ های مختلف حساب بشه و ببینیم هر وبلاگ با چند تا کلمه مستقل از هم نوشته شده (: البته معلومه که «می‌رود، می رود، میرود» چهار کلمه جدا شمرده شدن… برای حل نسبی این مشکل توی اون خط که رجکس یکسری کاراکتر اضافی رو حذف می کنه، نیم‌فاصله رو هم اضافه می‌کنم و خروجی اینطوری می‌شه:

jadi@jeducation:~/Downloads/weblog_word_usage$ ./count_all_words.py | wc -l
23405


جالب نیست؟‌ برای چهار هزار و چهارصد کلمه، من گاهی از نیم فاصله استفاده کردم و گاهی نکردم. حالا مهم نیست (:‌ قدم بعدی سورت کردن است. با جواهر گنو و نشون دادن بیست تا بالایی:

jadi@jeducation:~/Downloads/weblog_word_usage$ ./count_all_words.py | sort -n -r | head -20

خروجی رو براتون نمی‌ذارم چون فقط یکسری حرف ربط بی ربط است.

پروژه جانبی بسیار مهم برای زبان فارسی: همه زبان‌ها یک فایل دارن به اسم نمی دونم چی (کسی می‌دونست لطفا بگه) که توش کلمات «بی ربط» اون زبان نوشته شدن. به اصطلاح همون am و is و are یا امثال «است» و «شد» و «و» و «یا» و … که در اینجور جاها کاربرد داره (می شه اون کلمات رو از فهرست این کلمات که برنامه بهمون داده حذف کرد تا کلمات اختصاصی من به دست بیاد و نه چیزهایی که برای جمله ساختن همه استفاده می کنن). آیا چنین فایلی داریم برای زبان فارسی؟ تو خارجی‌ها اسمش چیه؟

حالا که بیست تا اولی به درد نخوردن، به جاش کل لیست رو می‌برمش توی لیبره آفیس و مثل همیشه نمودار بیشترین کلمات (غیر حرف ربط و استاندارد مثل است و باشد و شد و اینها) رو می‌کشیم تا با یک نمودار قشنگ کار رو تموم کرده باشیم:

پ.ن. این مجموعه مهمل ، یک قسمت دیگه هم داره (: کمی دیرتر ولی.